Mon. Not. R. Astron. Soc. 000, [HST] (2005) Printed 28 July 2008 (MN I^TeX style file v2.2) 

The turbulent destruction of clouds - I. A k-e treatment of 
turbulence in 2D models of adiabatic shock-cloud 
interactions 
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ABSTRACT 

The interaction of a shock with a cloud has been extensively studied in the litera- 
ture, where the effects of magnetic fields, radiative cooling and thermal conduction 
have been considered. In many cases, the formation of fully developed turbulence has 
been prevented by the artificial viscosity inherent in hydrodynamical simulations. This 
problem is particularly severe in some recent simulations designed to investigate the 
interaction of a flow with multiple clouds, where the resolution of individual clouds is 
necessarily poor. Furthermore, the shocked flow interacting with the cloud has been 
assumed to be completely uniform in all previous single-cloud studies. In reality, the 
flow behind the shock is also likely to be turbulent, with non-uniform density, pressure 
and velocity structure created as the shock sweeps over inhomogenities upstream of 
the cloud (as seen in recent multiple cloud simulations). To address these twin is- 
sues we use a sub-grid compressible k-e turbulence model to estimate the properties 
of the turbulence generated in shock-cloud interactions and the resulting increase in 
the transport coefficients that the turbulence brings. A detailed comparison with the 
output from an inviscid hydrodynamical code puts these new results into context. 

Despite the above concerns, we find that cloud destruction in inviscid and k-e 
models occurs at roughly the same speed when the post-shock flow is smooth and when 
the density contrast between the cloud and inter-cloud medium, x ^ 100- However, 
there are increasing and significant differences as x increases. The k-e models also 
demonstrate better convergence in resolution tests than inviscid models, a feature 
which is particularly useful for multiple-cloud simulations. 

Clouds which are over-run by a highly turbulent post-shock environment are de- 
stroyed significantly quicker as they are subject to strong "buffeting" by the flow. 
The decreased lifetime and faster acceleration of the cloud material to the speed of 
the ambient flow leads to a reduction in the total amount of circulation (vorticity) 
generated in the interaction, so that the amount of vorticity may be self-limiting. Ad- 
ditional calculations with an inviscid code where the post-shock flow is given random, 
grid-scale, motions confirms the more rapid destruction of the cloud. 

Our results clearly show that turbulence plays an important role in shock-cloud in- 
teractions, and that environmental turbulence adds a new dimension to the parameter 
space which has hitherto been studied. 

Key words: hydrodynamics - ISM: clouds - ISM: kinematics and dynamics - shock 
waves - supernova remnants - turbulence 
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1 INTRODUCTION 

Circumstellar, interstellar, and intergalactic environments 
are inhomogeneous, with clouds of various densities and tem- 
peratures embedded in a hotter, more tenuous, substrate. 
This substrate is often turbulent, in part due to the interac- 
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tion of shocks, shells, winds and jets with these clouds. The 
nature of such interactions is interesting, because the evolu- 
tion and morphology of large-scale flows can ultimately be 
regulated by objects of much smaller size. 

In the interstellar medium (ISM), for instance, the 
interaction of supernova shock waves with interstellar 
clouds creates a continuous interchange of mass and en- 



ergy between various th ermal phases (|Cox k. Sm ith 1974; 
iMcKee fc Ostrikerl 119771 ). Several outcomes are possible. 
The shocks may destroy the clouds, and mix their mate- 
rial into their surroundings. Alternatively, the shocks may 
trigger the collapse of the clouds and the formation of new 
stars, the reby removing material (a t least temporarily) from 
the ISM (lElmegreen fcLadal [19771 ). The shocks themselves 
will slow and material behind them will cool to form thin 
dense shells. These shells may then fragment and form new 
clouds. Clouds which survive the passage of the shell sub- 
sequently find themselves either inside a hot, low density 
bubble (if the bubble is energy-conserving), or exposed to a 
fierce, high Mach number wind (if the shell is momentum- 
driven) . An important process is the generation of vorticity 
as the diffuse medium flows around the clouds. Understand- 
ing the interaction between shocks/shells/ winds/jets and in- 
terstellar clouds is therefore a key step in studies of the 
str ucture and evolution of the ISM (see the recent re views 
bv lElmegreen k ScalollJOOJ : IScalo fc Elmegreenll2004l ). 

A similar interaction occurs in massive, early- type, stel- 
lar systems where each star blows a powerful, clumpy, wind 
which collides with the other. The impact of the clumps 
causes the wind- wind collision region to become highly tur- 
bulent (|Pittardll2007al l , with implications for particle accel- 
eration, the timescales for equilibrium ionization and elec- 
tron heating, and the physical mixing of the winds. 

In this work we consider the adiabatic interation of a 
shock with a cloud. An extensive literature of analytical and 
numerical in vestigations of shock-cl oud interactions now ex- 
ists (see, e.g. jNakamura et al.ll2006l . and references therein). 
The effects of ordered magnetic fields, radiative cooling, 
and thermal conduction have all been considered, but not 
simultaneously until recently (j Orlando et al.i i2QQ8i) . High- 
power laser expe riments of shock-cloud interactions (e.g., 
iKlein et aUl2003l ) have complemented this literature. While 
thermal conduction acts to prevent hydrodynamic instabili- 
ties, radiative cooling enhances them. The effect of ordered 
magnetic fields is more complicated: instabilities are pre- 
vented in cases where a magnetic field provides a high ten- 
sion at the surface of the cloud, but the effects depend on 
the strength, orientation and scale of the magnetic field. In- 
stabilities and vortical motions may be prevented in some 
directions, but not necessarily in others. Interestingly, in 
three-dimensional simulations a magnet ic field may actually 
enhance the fragmentation of a cloud (jGregori et al.lll999l : 
IShin. Stone fc Snyder 2008), though the actual mixing of 
the cloud and ambient medium remains hindered. Even if 
the magnetic fields are not strong enough to directly affect 
the dynamics, they can significantly reduce the effects of 
thermal conduction in directions normal t o the field line s, 
thus allowing instabilities to develop (Orla ndo et al.ll2008l V 

Given the fact that in shock-cloud interactions the 
Reynolds number is typically high, and the above predis- 
position for instabilities to develop, the fluid velocity field 
around the cloud is expected to vary significantly and irreg- 



ularly in both position and time. This "turbulence" trans- 
ports and mixes the cloud material much more effectively 
than a comparable laminar flow. Turbulence is also effective 
at "mixing" the momentum of a fluid (i.e. accelerating ma- 
terial ripped off the cloud to the ambient flow speed) , and at 
transferring heat. In turbulent flows of high Reynolds num- 
ber there is a separation of scales - large-scale motions will 
be strongly influenced by the geometry of the cloud, and 
control the transport and mixing of cloud material, while 
the behaviour of small-scale motions is determined almost 
entirely by the rate at which they receive energy from the 
large scales, and by the viscosity. Hence these small-scale 
motions have a universal character, independent of the flow 
geometry. 

While much insight has been gained from previous nu- 
merical investigations of shock-cloud interactions, the arti- 
ficial viscosity inherent in all such simulations has the po- 
tential to prevent the formation of fully developed turbu- 
lence, to limit the turbulent mixing of cloud material into 
the surrounding medium, and to hinder the destruction of 
the cloud. Furthermore, all previous simulations are highly 
idealized in the sense that the flow behind the shock is as- 
sumed to be perfectly smooth and uniform (i.e. laminar). 
In reality, random inhomogenities in the ambient medium 
upstream of the cloud will deform the shock, and will cause 
velocity, density and pressure structures to develop in the 
post-shock flow. The post-shock flow will then also be "tur- 
bulent" , and it is expected that the destruction of the cloud 
will be more rapid in such conditions. 

In this work we address whether fully developed turbu- 
lence has been prevented in all previous numerical works. For 
high Reynolds number flows, the only tractable method is a 
statistical approach i.e. to describe the turbulent flow, not in 
terms of a velocity field u(x, t), but in terms of some statis- 
tics. A model based on such statistics can lead to a tractable 
set of equations, because statistical fields vary smoothly (if 
at all) in position and time. To do this we use a sub-grid 
turbulent viscosity model, where an attempt is made to cal- 
culate the properties of the turbulence and the resulting 
increase in the transport coefficients. The most widely used 
is the so-called k-e model, where the properties of the tur- 
bulence are described by two variables, the turbulent energy 
per unit mass, /c, and the turbulent dissipation rate per unit 
mass, e. The addition of viscous and diffusive terms in the 
fluid equations simulates the turbulent mixing of the cloud 
and intercloud medium due to shear instabilities. In this 
way the subgrid turbulence model emulates a high Reynolds 
number flow. Incorporating a k-e model into shock-cloud 
simulations should produce more realistic results than those 
from inviscid codes, where the viscosity is purely numerical 
and the size of shear instabilities is determined by the res- 
olution of the numerical grid. A turbulent viscosity model 
differs from simply adding physical viscosity to the grid, be- 
cause the turbulent viscosity is largest in shear layers and es- 
sentially vanishes in regions with little shear, whereas models 
with grid viscosity have the same viscosity everywhere. 

The structure of this paper is as follows. The key 
physics of a shock-cloud interaction is reviewed in Section O 
Section [3] introduces the investigation and the numerical 
method used. The results are presented in Section 2] where 
the effects of turbulence on the cloud evolution are de- 
scribed. A discussion of the relevance of our results to shock- 
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cloud and wind-cloud observations is given in Section [5] Sec- 
tion [6] summarizes the conclusions of this work, and possible 
future work is noted in Section [71 



sity p 

and velocity U{ 



2.6 X 10"^^ g cm"^, temperature T : 



: 3.2 X 10^ K, 



0.75i;b 



1.1 X lO'' cms"'. The pa- 



rameters noted above are typical of small interstellar clouds 
dMcKee fc Ostriker 1977| ). 



2 THE INTERACTION OF A SHOCK WITH A 
CLOUD 

The interaction of a shock with a cloud is a highly non-linear 
and complex phenomenon which can be vastly simplified if 
some assumptions are made. In this work we assume that 
the magnetic field is too weak to be dynamically important 
(though it must be strong enough to reduce the thermal 
conductivity and effective mean- free-path) , and that the in- 
teraction is adiabatic. In the interstellar medium, the typical 
magnetic field is ^ 5 /iG, and the magnetic pressure is typi- 
cally a few times higher than the thermal pressure fe.g. JCoxl 
I2OO5I ) . However, if the shock driven into the cloud is strong 
and adiabatic the postshock thermal pressure becomes much 
higher than the postshock magnetic pressure, justifying our 
assumption that the magnetic field is too weak to be dynami- 
cally important. We also assume that the magnetic field does 
not prevent the full mixing of initially disparate phases (i.e. 
that turbulence drives the reconnection needed to allow this 
mixing) . 

The effects of thermal conduction are also ig- 
nored in this work. The efficiency of thermal conduc- 
tion in magnetized turbulent plasm as remains highly 
uncertain (see, e.g., IPittardI l2007bh . X-ray observa- 
tions of the hot intraclu ster medium ([Ettori fc FabianI 
l200d: IVikhlinin et al.l l200lh and theoretical considerations 
(iNaravan fc MedvedevI l2Q0ll: lAsai. Fukuda fc Matsumotol 
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I2OO4I : IChandran fc Maronll2004[ ) suggest that the coefficient 
of heat conduction is at least five times lower than the 
Spitzer value in the presence of tangled magnetic fields. 
Furthermore, if turbulent resistivity can reconnect field 
lines quickly enough, turbulent heat tr ansport niay be 
more efficient than thermal conduction (|Cho et al.l l2003l : 
IChandran fc Maronll20o3 : lLaza rian"2QQ6-V 

For the interaction to be adiabatic, radiative cooling 
must be unimportant. This is often the case for clouds in hot, 
low density environments (e.g., in planetary nebulae, in bub- 
bles blown around individual or groups of massive stars, and 
in starburst and superwind environments). The behaviour 
of the cloud is more likely to be adiabatic if the cloud is 
small. Since the assumption of adiabaticity preserves the 
scale-free nature of the simulations, the interaction is then 
determined by whether there is a sufficient rate of collisions 
within the plasma to give it fluid properties, by the domi- 
nant mechanism for the damping of hydromagnetic waves, 
by the Reynolds number, and by the way the turbulence is 
driven. Each of these issues is discussed below. In order that 
their importance can readily be determined, two different 
examples are considered. In the first scenario the cloud is 
ionized, and has a radius re — 1 pc, density ric = 0.4 cm""^, 
and temperature Tc = 8000 K. In the second scenario we 
consider a neutral cloud with Tc = 0.46 pc, ric — 30 cm""^, 
and Tc = 100 K. In both cases we imagine that the clouds 
are in approximate pressure equilibrium with a surrounding 
medium with nic = 3 x 10~^ cm~^ and Tic = 10^ K, and are 
struck by a high speed Mach 10 shock with velocity Vh — 
1.52 X 10^ cms~^. The shocked intercloud gas has a den- 



2.1 Collisional or collisionless? 

The first consideration is whether the mean-free-path is 
short enough that a collisional treatment can be adopted 
(i.e. whether the plasma has fluid- like properties). If this is 
the case, complications due to collisionless wave-particle in- 
teractions and their associated effe cts such as Landau damp- 
ing can be ignored (|Parkerlll979l V The damping of hydro- 
magnetic waves in a collisionless plasma is stronger than in 
a collisional plasma, since collisions serve to suppress the 
wave-particle interactions which damp the wave so strongly 
in the collisionless case. 

The Coulomb collision cross-section for electron or ion 
scattering is a ?^ 10""^^ /TJ^ cm^, where Tev is the temper- 
ature in electron volts. For neutral material the gas atomic 
cross-section is '^ 10~^^ cm^. Hence the mean- free-path 
within the ionized and neutral clouds is A = l/na ^ 10^^ cm 
and ?^ 10^^ cm respectively. These values are much smaller 
than the cloud radii, so a fluid-like treatment is appropriate. 
On the other hand, the mean- free-path in the post-shock in- 
tercloud medium is A ~ 6 x 10^° cm, which is significantly 
larger than the radius of the clouds. However, even a small 
B-field will significantly reduce the effective mean-free-path 
(the mean-free-path is then of order the gyroradius, which 
for B — ?> /iG and T = 10^ K is rg ~ 10^ cm for protons), 
so the shocked intercloud gas can also be considered to be 
collisional. 



2.2 Reynolds number, eddies and instabilities 

The Reynolds number of flow past a cloud is Re = urdv^ 
where u is the average flow velocity past the cloud, Tc is the 
radius of the cloud, and v is the kinematic viscosity. For a 
fully ionized, non-magnetic gas of density p and temperature 



2.21 X 10" 



_ ^5/2^1/2 
Z^p In A 



(1) 



where A and Z are the atomic weight and charg e of the 
posit ive ions, and In A is the coulomb logarithm (|Spitzerl 
ll956l V In a magnetized plasma, the kinematic viscosity is of 
the order of the mean- free-path (i.e. the particle gyroradius) 
times the typical thermal velocity. 

The characteristic Reynolds number of the interaction 
is higher in the cloud material than the surrounding envi- 
ronment, since the cloud is considerably cooler. The shock 
driven into the cloud has a speed v — v^/ ^/x^ where x is the 
density contrast of the cloud with respect to its surround- 
ings. For the ionized cloud considered in Section [21 the post- 
shock temperature is 2.6 x 10^ K and z/ ?^ 5 x 10^° cm^ s""*^. 
Hence Re ~ 7 x 10^. If the cloud is magnetized with a post- 
shock field B r^ 1^ /iG, the kinematic viscosity for protons 
is zy ~ 10^^ cm^ s~^ and the Reynolds number is propor- 
tionally higher. In the neutral cloud the Reynolds number 
calculated from the kinematic molecular viscosity is simi- 
larly high. The Reynolds number of the flow around a cloud 
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is also high. For a flow with B — ?> /xG and T — 10^ K, 
the kinematic viscosity v ^ 10^^ cm^ s~^. In all cases, the 
viscous stresses acting on the boundary layer which forms 
as the shock sweeps over the cloud are negligible, and a tur- 
bulent energy cascade ensues. 

The largest eddies have a length scale, /, comparable in 
size to the cloud, while the smallest eddies where the turbu- 
lent energy is dissipated have a length scale 77 ~ Re""^/^ /. 
Due to the nonlinear term u • Vu in the equation of mo- 
tion, large eddies, created by instabilities in the mean flow, 
are themselves subject to inert ial instabilities and rapidly 
"break-up" or evolve into yet smaller vortices. Energy is 
transferred into vortices of about one- half their size in a 
time comparable to their "turnover time" (t — l /u, where 
u is t he characteristic velocity of eddies of size /, IPavidsonl 
I2OO4I ). Provided that energy is constantly injected at large 
scales, small-scale eddies are superimposed on larger eddies. 
The timescale to set up a turbulent energy cascade is roughly 
twice the turnover time of the largest eddies (i.e. t — 21 /u). 
Since l/u ^ rc/i^ics, the setup time is t ^ 2rc/vh = tsc, 
where tsc is the timescale for the shock in the intercloud 
medium to sweep over the cloud. For dense clouds (x ^ 1)? 
this setup time is much shorter than the survival time of the 
cloud, which is a few times the "cloud crushing" timescale. 



_ X^^Ve 



Wb 



(2) 



for the cloud to be crushed by the init ial shock that is driven 
into it (JKlein. McKee fc Colellalll994l ). 

Simultaneously with the top-down energy transfer from 
large to small eddies, Kelvin-Helmholtz (KH) and Rayleigh- 
Taylor (RT) instabilities inject energy from the bottom-up, 
since the smallest scale disturbances grow fastest. The KH 
and RT growth times are 



tKH 



k\rc 



tRT ' 



(fcArc)i/2' 



(3) 



where k\ is th e wave-number of the perturbation 
(|Klein et al.l 1 19941 ) . The smallest scale of the instabilities 
is set by the scale at which the damping of hydromagnetic 
waves occurs, which was shown in Section [ZTI to be through 
particle collisions rather than wave-particle interactions. 

For the unmagnetized ionized cloud the minimum scale 
due to viscous damping is ?7vis '^ Re""^/^ Tc '^ 4 x 10~^ Vc- 
However, in ionized plasmas the thermal conductivity is even 
more effective at damping waves with a significant longitudi- 
nal component. The thermometric conductivity K — HiT /U ^ 
where U is the thermal energy per unit volum e and Hi is the 
thermal conduction coefficient (|Parkerlll979l l. Here we find 
that K ^ AQ v^ and thermal conduction prevents instabili- 
ties with a scale smaller than ?7tc '^ 6 x 10~^ Vc. In magne- 
tized ionized clouds the lengthscale at which instabilities are 
damped is smaller if the magnetic field is sufficiently weak 
to be dynamically unimportant. On the other hand, if the 
magnetic field is strong enough to be dynamically impor- 
tant, magnetic tension increases the minimum lengthscale 
at which instabilities occur. 

A disturbance of wavelength A = 10 ""^ Tc grows in a 
timescale of txH = 1.6 x 10~^ tec and tRx = 1.3 x 10~^ tec- 
With the ionized cloud parameters given above, x = 133, 
tec = 5.8 tsc, and we obtain Iku ~ 10""^ tsc and tRx ~ 
0.1 tsc- These timescales are much faster than the one to 



setup the turbulent energy cascade (which was shown above 
to be ^ tec), and the (ensemble- averaged) turbulent spec- 
trum will differ from the classical Kolmogorov spectrum be- 
cause it is driven by energy input at both large and small 
scales. 

Photoionization at the surface of neut ral clouds main- 
tains most of the C in the form C~^ (e.g., iHartquist et al.l 
"1998"). An upper limit to the ratio of the ionized to neutral 
mass density, p*/p, is therefore the fractional abundance by 
mass of carbon, which has a cosmic value of about 3 x 10~^ 
fe.g. JPopita fc Sutherlandll2003l ). lKulsrud k Pearcd (Il969l ) 
note that for hydromagnetic waves with A < A2 , the charged 
particles move as if the neutrals were absent, while for 
A > Ai, the entire medium moves. In contrast, when A2 < 
A < Ai, the neutrals and ions move independently and fric- 
tion is strong, and waves do not propagate. Ai corresponds 
to waves with angular frequencies comparable to the rate at 
which neutrals transfer momentum to ions, while A2 is the 
corresponding rate for momentum transfer from ions to neu- 
trals. Assuming that the magnetic field strength within our 
neutral cloud is ^ 3 /xG, we find that Ai ~ 4 x 10^^ cm and 
A2 ~ 1.3 X 10^^ cm. In this case fully developed turbulence 
is prevented, since Xi/rc > ?7vis when Re ^ 10^, although 
instabilities continue to be driven from the bottom up. Fully 
developed turbulence may be obtained when the magnetic 
field is weaker, since Ai and A2 are both proportional to B, 
or when the cloud is larger. 

In numerical simulations of shock-cloud interactions, 
the growth of KH and RT instabilities is closely related to 
the development of the slip surface around the cloud, and 
perturbations with wavelengt hs smaller than the thic kness of 
the shear layer are stabilized (|Nakamura et al.ll2006l V Thus, 
small scale instabilities have been artificially suppressed in 
previous work on hydrodynamic shock-cloud interactions. 
Our use of a /c — e turbulence model in this paper attempts 
to address this shortcoming. 



2.3 The turbulent boundary layer 

IHartquist fc DvsonI (|l988r ) argued that the turbulent bound- 
ary layer which forms around a cloud has a thickness of 
the order of rc/VRet, where Ret is an effective "turbulent" 
Reynolds number arising from the fact that the turbulence 
itself gives rise to an effe ctive viscosity. Since Ret ~ 10^ 
(jHartauist fc DvsonI Il988l ) , the thickness of the turbulent 
boundary layer is a few percent of the cloud radius. More 
detailed calculations and laboratory experiments show that 
the opening angle of a turbulent mixing layer in a mildly 
supersonic flow (such a s occurs behind a hig h Mach number 
shock) is of order 10° (|Cant6 fc Ragalll99lh - we find good 
agreement with this (see Section 14. 2p . Convergence tests re- 
veal that a minimum numerical resolution of about 120 cells 
per cloud radius is needed for convergence of various global 
quantities (see, e.g., Nakamura et al. 2006, although the nec- 
essary resolution may to some extent also depend on the 
numerical scheme). This is consistent with simulations of 
this resolution and higher beginning to resolve the turbu- 
lent boundary layer. 
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2.4 Cloud destruction and mixing 

The main stages in the destruction of a cloud by a shock and 
the subsequent mixing of its material into the surrounding 
flow are reviewed in Section 14.11 The most disruptive KH 
and RT instabilities are those with wavelengths of the or- 
der of the cloud radius. However, there are a number of 
ways in which the growth of KH and RT instabilities may 
be hindered or amplified. For instance, clouds with diffuse 
boundaries are less susceptible to KH instabilities and sur- 
vive longer (INakamura et al.ll2006[ ). Two dimensional MHD 
simulations found that the growth of KH and RT instabili- 
ties is strongly inhibited when there is a dynamically impor- 
tant magnetic field, due to the field providing an additional 
tension at the interface betwee n the cloud and the surround- 
ing flow ([Mac Low et al.ll 19941 ) . Somewhat surprisingly, fully 
three dimensional magneto- hydrodynamic (MHD) simula- 
tions of a wind-cloud interaction revealed that an ordered 
magnetic field can actually enhance hydrodynamic instabil- 
ities, as background field lines b ecome trapped in d eforma- 
tions in the surface of the cloud ([Gregori et al.lll999l . though 
these were low resolution calculations). Higher resolution 
shock- cloud simulations recently presented by IShin et al.l 
([20081) show that irrespective of the field geometry, and the 
morphology of the cloud fragments which are produced in 
the interaction, the rate of mixing is reduced compared to 
the non-magnetic case. However, sufficiently weak magnetic 
fields have no dynamical influence, and only offer a potential 
reduction in the collision mean- free-path. Strong thermal 
conduction supp resses hydrodynamic instabilities (see Sec - 
tion [2:21 and also lOrlando e t al. 2QQ5l: lMarcolini et al.ll2005l ). 
but the degree of t his effect is sensitive t o the orientation of 
any magnetic field (|Orlando et al.ll2008l V 

KH and RT instabilities are always stronger when radia- 
tive cooling is important, and the cloud breaks up into nu- 



merous dense, cold fragments JMellema, Kurk fc Rottgerind 
I2OO2I : [Fragile et al]l2004l . 120051 : IVan Loo et al.ll2007l ). In the 



corresponding simulations the fragments appear to survive 
for an appreciable time, but are poorly resolved, so the 
timescales corresponding to their further evolution are some- 
what uncertain. The way in which clouds are destroyed and 
mixed into their surroundings is sensitive also to the density 
contrast between the cloud and the surrounding medium, in 
that clouds with higher values of x suffer direct stripping 
of material from their surfaces by hydrodynamic ablation 
(JKlein et al.1 Il994 also compare Figs. [H [9l and [10] in this 
work) . 

Cloud destruction and the mixing of the cloud mate- 
rial into the surrounding flow are two distinct processes 
which were not always distinguished in previous work. In 
INakamura et al.l (|2006r ), the cloud destruction timescale is 
taken to be the time when the largest fragment drops be- 
low a certain fraction of the initial cloud mass, while the 
timescale for the mixing of former cloud material into the 
surrounding medium is estimated by comparing the inte- 
grated mass above a particular threshold density with the 
initial cloud mass. However, even weak magnetic fields can 
prevent the actual mixing of stripped material with the sur- 
rounding medium, and reconnection on small-scales is nec- 
essary if the plasmas are to mix fully. This may, of course, 
take some considerable time to achieve. 



3 THE NUMERICAL SETUP 

3.1 The numerical scheme 

The shock-cloud interaction is modelled by solving numer- 
ically the Euler equations of inviscid fluid flow, supple- 
mented by a sub-grid turbulent viscosity model as appro- 
priate. When the sub-grid model is included, the continuity, 
scalar, momentum, energy, turbulent energy and turbulent 
dissipation equations are respectively: 



| + v.(.„, 
4^ + v.,..„) 



0, 



^ + V-(puu) + VP 
dE 



V-r: 



dt 



+ V-[(^ + P)u-u-t] 



|S.v.(,„, 



0, 

7-1 
V • i/irVk) = Sk, 

V • (M.Ve) = S,. 



(4) 
(5) 
(6) 

{iitVT) = Se, (7) 

(8) 
(9) 



Here p is the mass density, u is the velocity, E is the total 
energy density 

P 



7^ + ^^"'' 



(10) 



P is the thermal pressure, k is the turbulent energy per unit 
mass, e is the turbulent dissipation rate per unit mass, and 
the turbulent diffusion coefficients are 



Mt 






IJ-^ = 



1.3' 



(11) 



where C^ = 0.09. a^: is an advected scalar used to distinguish 
between cloud and ambient material. 

The momentum equation source term, Sp, is zero in 
Cartesians. In cylindrical symmetry it is 



/ 


ff^ — ^^1 

.3r r^ . 


r 


y+ipk 


\ 


■ (12) 


V 


/ 




The k and e source terms are respectively 






Sk^Pt- pe, 




(13) 


and 






Se^^iClPt 


-C2pe), 








(14) 



where Ci = 1.4 and C2 = 1.94. 

The turbulent production term 



Pt = PT 



dui 
dxj 



dui du 
dxj dx 



fV-u(p/c + /iTV-u),(15) 



where the summation convention is assumed. The terms in- 
volving /iT are due to the rate of working of the turbulent 
stresses, and the terms involving V • u take into account 
volume changes on the turbulence. In cylindrical symmetry, 
the production term has to be complemented by an extra 
geometric term 

2^T^. (16) 
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The turbulent stress tensor, r, is 



^ MT 



duj 
dxj 



duj_ 
dxi 



-^ij(/iTV- u + p/c). 



(17) 



All computations were performed in 2D cylindrical sym- 
metry using an Eulerian adapative mesh refinement (AMR) 
hydrodynamic code, with a linear Godunov solver and piece- 
wise linear cell interpolation fsee lFalldll99lh . Although it is 
of lower order than the piecewise parabolic method (PPM), 
it performs well in multi-dimensional problems, as it is only 
partially operator split, and it performs better than PPM 
for problems wh ere there is rapid advectio n across the com- 
putational grid ([Runacres k, Owockil 120051 ) 

The entire computational domain is covered by the two 
coarsest grids, G^ and G^ . The solution at each position 
is calculated on all grids that exist there, and the differ- 
ence between these solutions is used to control refinement 
(note that refinement criteria based on the local gradients 
of only selected variables (e.g., de nsity) do not properly re- 
solve turbulent flow - see lapichin o et al.ll2008r ). Finer grids 
are dynamically added where they are needed, and removed 
where they are not. Each refinement level increases the res- 
olution in each of the spatial directions by a factor of 2, 
and the refinement is done on a cell-by-cell rather than on a 
patch basis. The time-step on grid G^ is Ato/2n, where Ato 
is the time-step on C^, in order to ensure Courant number 
matching at the boundaries between coarse and fine grids. 

The k-e model is designed to model the mean flow 
in fully developed, high Reynolds number, turbulence. It 
has been calibrated by comparing the computed growth of 
shear layers with exp eriments of high Reynolds number flows 
(iDash k WolJIlQsd ). Although the flow in our problem is 
somewhat more complicated than in these experiments, it 
is described by exactly the same equations. Of course, any 
turbulence model can only be approximately correct, but it 
should give more reliable results than an inviscid calculation. 

In the sub-grid model, turbulent energy is generated 
by the action of the turbulent viscosity on the mean flow 
and is converted to heat at the dissipation rate, e. Since 
the turbulent energy resides mainly in large eddies, while 
the dissipation occurs in the small ones, one can think of 
k and e as describing the large-scale and small-scale turbu- 
lence respectively. Since the aim of the sub-grid model is 
to mimic a three dimensional turbulent flow, the effects of 
the turbulence (such as enhanced transport coefficients) are 
treated correctly, even though the grid is cylindrically sym- 
metric. However, the turbulent motions that are resolved on 
the grid are actually vortex rings, and not eddies. F urther 
details of the model implementation can be found in iFalld 
(|l994l ). 



3.2 Initial and boundary conditions 

We consider a Mach 10 shock interacting with a cloud with 
a density contrast x of either 10, 10^ or 10^, and compute 
simulations with different spatial resolutions using either an 
inviscid code or one which includes a k-e turbulence model. 
The effect of different levels of turbulence in the post-shock 
gas is also explored. Table [1] summarizes the calculations 
performed. Most computations are for clouds with steep 
density profiles, but we have also examined the effect of a 



Table 1. Summary of the shock-cloud simulations performed. 
The resolution is the number of cells per cloud radius on the finest 
grid. Models with "sh" in their name were computed for clouds 
with a shallow density gradient. Model clrtb32 was computed 
using an inviscid code with grid-scale turbulence. 



Model 


X 


Pi 


resolution 


turbulence 


clno 


IQi 


10 


32,64,128 


no 


cllo 


10^ 


10 


32,64,128 


low 


clhi 


IQi 


10 


64,128 


high 


c2no 


10^ 


10 


32,64,128,256 


no 


c21o 


10^ 


10 


32,64,128,256 


low 


c2hi 


10^ 


10 


16,32,64,128,256 


high 


c3no 


10^ 


10 


16,32,64,128 


no 


c31o 


10^ 


10 


16,32,64,128 


low 


c3hi 


103 


10 


32,64,128 


high 


c2nosh 


10^ 


1 


64 


no 


c21osh 


102 


1 


64 


low 


c2hish 


102 


1 


64 


high 


c3nosh 


103 


1 


64 


no 


c31osh 


10^ 


1 


64 


low 


c3hish 


10^ 


1 


64 


high 


clrtb 


lOi 


10 


32 


grid-scale 



shallower density profile on the resulting evolution (see Sec- 
tion l4.4p . A calculation using an inviscid code with grid-scale 
post-shock turbulence is also presented (see Section |43)) . 

The calculations are computed on an r — z cylindrically 
symmetric grid, with a domain of ^ r ^ 24, —94 ^ 2; ^ 6 
when X = 10, ^ r ^ 24, -120 ^ z ^ 6 when x = 10^ and 
^ r ^ 48, -594 ^ z ^ 6 when x = 10^. This ensures that 
the cloud is well dispersed and mixed into the post-shock 
flow before the shock reaches the edge of the numerical grid. 
In this way the global quantities detailed in Section 13.41 are 
accurately computed. All calculations are for an ideal gas 
with 7 = 5/3, and are scaled so that the fluid variables have 
values reasonably close to unity. 

Several additional parameters must be specified when 
using a turbulence model. One of these is the maximum 
eddy size, which here is set equal to the cloud radius. An- 
other choice concerns the initial level of turbulence in the 
gas, and two extreme cases are considered. In the first case 
(hereafter identified by "low /c-e" , or "lo" in the model name) 
the postshock gas initially has an extremely low level of tur- 
bulence, with a ratio of turbulent energy density to thermal 
energy density etb/eth ^ 10~^. In the second case (identi- 
fied by "high /c-e", or "hi" in the model name) this ratio is 
0.13 (values much higher than this cause the shock to ac- 
celerate too much and the interaction does not occur at the 
intended Mach number). High levels of post-shock turbu- 
lence may arise if the shock is propagating through an inho- 
mogeneous medium (e.g., when there are density variations 
further upstream). As we shall see, this ratio is similar to 
the turbulent energy fraction attained by the cloud material 
after the shock encounter (see Section l4.6.7p . so there is a 
degree of self-consistency in these models. In both cases, ini- 
tially etb/eth = 0.04 within the cloud (this is set low enough 
to not affect the dynamics - cf. the initial development of the 
interaction in models c3nol28 and c31ol28), and the cloud 
and inter cloud medium are set in pressure equilibrium. 
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Figure 1. Comparison of cloud density profiles obtained using 
Eq. [18] (pi =1 and 10) and Eq. 1 in iNakamura et al.l (|2006l ) 
(n = 2 and 24) with x = 10^. 



3.3 Cloud density profile 

Clouds in the ISM do not have infinitely sha rp edges (see, 
e.g., the discussion in INakamura et al.l 120061 ). The density 
profile adopted in this work is 



p{r) = pamhl'ip + {I -ip)r]], 



where 

1 

^=2 



1 + 



a+ 1 



(18) 



(19) 



a = exp {min[20.0,pi((r/rc)^ — 1)]}, and r is the distance 
from the centre of the cloud of radius re. '0 is adjusted to 
obtain a specific density contrast for the centre of the cloud 
with respect to the ambient medium, x — Pmax/pamb- The 
parameter pi controls the steepness of the profile at the edge 
of the cloud. Eq. [18] tends to give a flatter density profile 
within the centre of the cloud, and a steeper profile as the 
cloud merges into th e ambient mediu n i, tha n profiles ob- 
tained using Eq. 1 in INakamura et aU (|2006h . For pi > 5, 
^ ~ X- Clouds with reasonably sharp ed ges are obtained 
with p i = 10 (similar to those from Eq. 1 in INakamura et al.l 
((20061) with n = 24) , while pi = 1 p roduces a much more ex - 
tended cloud which is closer to the INakamura et al.l (|2006|) 
profile with n = 2 (see Fig. [T}. 



3.4 Global quantities 

The cloud ev olution is studied through various integ rated 
quantities fsee lKlein et al.lll994l : [Nakamura et al.ll2006[ i. Av- 
eraged quantities (/), are constructed by 



(/) 



npf dV, 



1 
where the mass identified as being part of the cloud is 



(20) 



mi3 ■ 



j ' 



Kp dV. 



(21) 



K is an advected scalar, which has an initial value of 
p/ixpamh) for cells within a distance of 2.25rc from the cen- 
tre of the cloud, and a value of zero at greater distances. 



Hence, k, = 1 in the centre of the cloud, and declines out- 
wards. The above integrations are performed only over cells 
in which k is at least as great as the threshold value, f3. 
Setting /3 = 0.5 probes only the densest parts of the cloud 
and its fragments, while setting (3 = 2/x probes the whole 
cloud including its low density envelope, and regions where 
only a small percentage of cloud material is mixed into the 
ambient medium. 

To measure the shape of the cloud, effective radii normal 
to and along the axis of symmetry are defined respectively 



"(jC'O 



1/2 



c=[^{{z^)-{zf)] 



1/2 



(22) 



In inviscid calculations a measure of the turbulence of the 
cloud is obtained from the velocity dispersions in the radial 
and axial directions, defined respectively as 

6v. = {v',Y'\ 5v. = ^(^yl)_^y^)y\ (23) 

The mean density is defined as 

W^f. (24) 

where Vjs is the volume of the region with n^ [3. 

All quantities computed with /3 = 0.5 are identified with 
the subscript "core" (e.g., acore), while those computed with 
(3 = 2/x are given the subscript "cloud" (e.g., adoud)- 

3.5 Timescales 

Several timescales are obtained from the simulations. The 
characteristic radial expansion timescale, tm, is defined as 
the time at which the cloud radius normal to the axis of 
symmetry, a, has increased to 90 per cent of its maximum 
value. The time for the cloud velocity relative to that of 
the postshock ambient flow to decrease by a factor of 1/e is 
defined as the "drag time", tdrag- The "mixing time", tmix, 
is defined as the time when the mass of the core of the cloud, 
^^core, reaches half of its initial value. The zero-point of all 
the time measurements quoted in this work occurs when the 
intercloud shock is level with the centre of the cloud. 



3.6 Convergence tests 

It is important to demonstrate that the calculations are per- 
formed at spatial resolutions that are high enough to resolve 
key features of the interaction. Increasing the resolution in 
inviscid calculations leads to smaller scales of instabilities. 
Quantities which are sensitive to these small scales (such 
as the mixing rate between cloud and ambient gas) may 
not be converged, while quantities which are insensitive to 
gas motions at small scales (e.g., the shape of the cloud) 
are more likely to show convergence. P revious studies (e.g., 
iKlein et alll 19941 : INakamura et al.ll2006l ) have indicated that 
about 100 cells per cloud radius are needed for convergence 
of the simulations. Here we carry out a similar study for 
calculations which use a subgrid turbulence model. 

Fig. [2] shows the evolution of the core mass and mean 
cloud velocity as a function of spatial resolution for both 
inviscid and k-e calculations. These parameters are a good 
test of the convergence, since both the rate of mixing and the 
momentum transfer between the cloud and the ambient flow 
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Table 2. Dependence of the global cloud and core properties on the level of turbulence and the density contrast of the cloud (see Table [T}. 
In each case the cloud was hit by a Mach 10 shock. The time-dependent quantities are evaluated at t = tmix rather than t = tm (c.f. 
iNakamura et al.ll20Q6l ). since adoud continues to rise in some simulations. Values in parentheses are obtained from integrations over the 
"core" mass rather than the "cloud" mass. Model names containing "sh" were computed using a shallow density profile {pi = 1). Model 
clrt32 was computed using an inviscid code and grid-scale post-shock turbulence. 



Model 



^drag/^c 



t^mix/^cc tm/tc 



a/rc 



c/Tc 



c/a 



{P)/Pn 



{V7.)/Vi, 



clnol28 


1.03 


(1.10) 


(6.82) 


3.82 


1.74 


1.72) 


2.82 


(2.92) 


1.62 


1.70) 


0.65 (0.87) 


0.634 


0.685) 


cllol28 


1.04 


(1.11) 


(5.96) 


3.57 


1.66 


1.74) 


2.20 


(2.40) 


1.33 


1.38) 


0.68 (1.03) 


0.579 


0.629) 


clhil28 


0.72 


(0.86) 


(4.37) 


3.52 


1.88 


1.90) 


1.36 


(1.30) 


0.72 


;o.68) 


0.56 (0.66) 


0.621 


0.651) 


c2nol28 


3.06 


(3.10) 


(5.05) 


4.53 


4.36 


4.21) 


5.07 


(3.58) 


1.16 


;o.85) 


0.073 (0.130) 


0.635 


0.606) 


c21ol28 


3.08 


(3.09) 


(4.86) 


4.82 


3.92 


3.35) 


4.12 


(3.05) 


1.05 


;o.9i) 


0.075 (0.142) 


0.593 


0.554) 


c2hil28 


2.47 


(2.65) 


(5.73) 


- 


4.05 


4.27) 


6.24 


(1.66) 


1.54 


0.39) 


0.058 (0.078) 


0.606 


0.590) 


c3nol28 


6.58 


(9.48) 


(8.59) 


10.23 


2.60 


1.22) 


75.9 


(12.6) 


29.2 


10.3) 


0.0081 (0.0401) 


0.394 


0.189) 


c31ol28 


6.84 


(7.15) 


(7.83) 


- 


4.01 


2.82) 


49.1 


(6.88) 


12.2 


2.44) 


0.0078 (0.0248) 


0.430 


0.271) 


c3hil28 


4.58 


(4.93) 


(5.95) 


7.25 


6.84 


7.79) 


38.2 


(5.57) 


5.58 


;o.7i) 


0.0064 (0.0122) 


0.584 


0.578) 


c2nosh64 


4.05 


(5.51) 


(7.92) 


8.60 


3.09 


2.75) 


11.7 


(2.33) 


3.80 


0.85) 


0.104 (0.392) 


0.455 


0.402) 


c21osh64 


4.00 


(4.33) 


(7.37) 


9.50 


2.89 


2.21) 


10.8 


(2.02) 


3.74 


0.91) 


0.103 (0.446) 


0.441 


0.377) 


c2hish64 


3.08 


(3.85) 


(4.87) 


- 


2.49 


1.46) 


6.50 


(0.98) 


2.61 


0.67) 


0.114 (0.476) 


0.379 


0.313) 


c3nosh64 


6.72 


(9.84) 


(7.86) 


- 


5.16 


2.42) 


48.6 


(8.64) 


9.41 


3.57) 


0.0124 (0.145) 


0.370 


0.216) 


c31osh64 


7.56 ( 


14.44) 


(8.99) 


5.16 


3.89 


1.05) 


64.2 


(6.12) 


16.5 


5.86) 


0.0121 (0.085) 


0.302 


0.150) 


c3hish64 


6.55 


(8.36) 


(4.95) 


11.82 


2.63 


2.13) 


22.8 


(2.94) 


8.67 


1.38) 


0.0241 (0.144) 


0.204 


0.155) 


clrt32 


0.66 


(0.72) 


(5.84) 


- 


2.21 


2.58) 


1.29 


(0.90) 


0.58 


0.35) 


0.61 (0.80) 


0.614 


0.666) 



c3no32 - 
c3no64 ■ 
c3no128 ■ 





Figure 2. Convergence tests for a shock-cloud interaction with x — 10^ ^ii^d pi = 10 using an inviscid (panels a and b) and k-e (panels c 
and d) code. The time evolution of the core mass (panels a and c) and the mean velocity of the cloud (panels b and d) are shown. Note 
the much tighter correlation of mcore and ("^^z, cloud) from the k-e code. 
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Figure 3. Relative error versus spatial resolution (number of cells per cloud radius on the finest grid) for a number of global quantities 
measured from a shock-cloud interaction with x = 10^ cind pi = 10 at t = 6 tec (left panels) and 8tcc (right panels). The top panels 
show results from the inviscid code, while the bottom panels are the results from the k-e code with "low /c-e" conditions (models c31o). 



are sensitive to small scale instabilities. It is therefore not 
too surprising to see that these quantities are poorly con- 
verged in the lo w resolution inviscid calculations (see also 
IShin et al.ll2008l V In contrast, the subgrid turbulence model 
leads to results which are much less dependent on the spatial 
resolution. This is also demonstrated in Fig. [3] which shows, 
for a number of parameters, the relative error defined as the 
fractional difference between the value measured at resolu- 
tion N and the value at the finest resolution, /: 



AQn 



\QN-Qf\ 
\Qf\ ■ 



(25) 



The convergence is much better for the k-e calculations (see 
Fig. [3]), leading to much less variation with resolution as seen 
in Fig. [21 Both inviscid and k-e calculations demonstrate 
better convergence at x = 10^ (not shown). 

Figs. [2] and [3] support claims from previous studies that 
of order 100 cells per cloud radius are needed for conver- 
gence. However, the neglect of the effect of turbulent eddies 
on the mean flow means that it is not clear that inviscid 
simulations, particularly at high values of x? are actually 
converging to a physically realistic solution. 



4 RESULTS 

4.1 Stages 

The main stages of the interaction of a shock of velocity Vh 
with a uniform cloud of density contrast x, in the adiabatic, 
un-magnetized, non-conducting case are: (1) an "initial tran- 
sient stage", where the incident shock propagates into the 
cloud with velocity Vs = Vh/x'^^^^ ^^^ ^ bow shock or bow 
wave propagates upstream into the ambient medium; (2) the 
"compression stage" , where the cloud is compressed mainly 
in the ^-direction by the transmitted shock and by a shock 
driven into the back of the cloud; (3) the "expansion stage", 
where the highly-pressured cloud expands downstream and 
laterally; and (4) the "destruction stage", where the cloud 
is destroyed and its material mixed into the surrounding 
flow. In other situations, for example when there is efficient 
cooling, the evolution can be significantly different (see Sec- 
tion. 12. 4[) . In all work to date the cloud is destroyed by the 
shock. The addition of gravitational forces is likely to be 
needed if the cloud is to survive. 



4.2 Cloud morphology and turbulence 

In Fig. m snapshots of the density distribution at different 
times are shown for an inviscid calculation with 128 grid 
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Figure 4. Snapshots of the density distribution from an inviscid calculation of a Mach 10 adiabatic shock hitting a cloud with a density 
contrast of 10^ with respect to the ambient medium and a density profile with pi = 10 (model c3nol28). The resolution is 128 cells per 
cloud radius. The evolution proceeds left to right and top to bottom with t = 0.0, 0.1, 0.3, 0.49, 0.87, 1.35, 1.83, 2.31, 2.79, and 3.75 tec- 
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Figure 5. As Fig. 0] but for a k-e calculation with low initial postshock turbulence (model c31ol28). The times of the snapshots are 
t = 0.1, 0.49, 0.87, 1.83, and 3.75 tec- 




Us 



Figure 6. As Fig. |4] but for a k-e calculation with high initial postshock turbulence (model c3hil28). The times of the snapshots are 
t = 0.1, 0.49, 0.87, 1.83, and 3.75 tec- 



12 J. M. Pittard, S. A. E. G. Falle, T. W. Hartquist and J. E. Dyson 



c3nol28 



C3I0128 



c3hil28 



4.0 6.0 



Figure 7. Comparison between the inviscid (model c3nol28, left panel) and k-e calculations with low and high initial postshock turbulence 
(models c31ol28 and c3hil28, middle and right panels, respectively) at t = 5.66 tec- 



cells per cloud radius (model c3nol28). The evolution of 
the cloud broadly follows the stages outlined above. The 
first two stages last until t ~ tec (i.e. the top 5 panels in 
Fig. |4|. The expansion of the cloud in stage 3 is supersonic 
with respect to the sound speed within the cloud, and a low 
density interior surrounded by a higher density shell forms 
(see the snapshot at t = 1.83 tec). The high density shell 
then collapses in on itself (see the snapshot at t = 2.31 tec)- 
Material is continuously ablated off the surface of the cloud 
by the fast-flowing surroundings, and a turbulent wake with 
prominent RT and KH instabilities forms. Fig. [71 shows that 
at later times the mass-loss from the cloud resembles a single 
tail-like structure (this is in contrast to models with lower 
values of X - see Section |43)) . 

Fig. [5] shows that the initial interaction of the shock 
with the cloud in the low k-e case (model c31ol28) is similar 



to the inviscid case, since the initial post-shock turbulence 
is low and k and e are small. However, the viscosity intro- 
duced by the sub-grid turbulence model prevents the sub- 
sequent development of the resolution-dependent RT and 
KH instabilities seen in Fig. [4] Instead, the loss of material 
from the cloud occurs much more smoothly. This is exactly 
as expected given that the purpose of the k-e model is to 
approximate the time- averaged flow. 

Simulations with a high initial level of post-shock turbu- 
lence are different again, as seen in Fig. [6] where the results 
of model c3hil28 are displayed. The transport/diffusion co- 
efficients are considerably higher in this simulation, and 
this leads to a faster rate of ablation from the cloud, and 
ultimately its more rapid destruction. The high level of 
upstream turbulence also smooths/broadens the bowshock 
ahead of the cloud and the tailshock which forms down- 
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Figure 8. Snapshots of the turbulent energy per unit mass, k, from the k-e calculations with x = 10^ cind pi = 10. The three left-most 
panels show the evolution with low initial postshock turbulence (model c31ol28) proceeding left to right with t = 0.1,0.3, and 0.87 tec- 
The two right-most panels show the evolution with high initial postshock turbulence (model c3hil28) proceeding left to right with t = 0.1 
and 0.87 tec The white regions in the middle panel are artifacts of the plotting routine. 



stream, so that their time- averaged positions are repre- 
sented. 

Fig. [71 compares the morphology of the clouds in these 
3 simulations at t = 5.66 tec The global features of mod- 
els c3nol28 and c31ol28 are reasonably comparable at this 
time, but the cloud in model c3hil28 is clearly at a more 
advanced stage of destruction (cf. Tableland Fig. I17|) . In- 
terestingly, material stripped from the cloud lies off-axis in 
model c3hil28. This develops from the off-axis density peak 
seen at t = 1.83 and 3.75 tec shown in Fig. [6l though the 
initial disturbance occurs at even earlier times. Clearly the 
turbulence in this model affects the properties of the shocks 
driven into the cloud and the slip surface that forms around 
it, with small differences at early stages being amplified dur- 
ing the subsequent evolution. 

The development of the sub-grid turbulence in model 
c31ol28 is shown in the three left-most panels of Fig. [H] 
where the turbulent energy per unit mass, k, is displayed. 
k is created in regions of high shear, particularly in a thin 
turbulent boundary layer along the slip surface around the 
cloud (see the left-most panel of Fig. (8]). The turbulent inten- 
sity quickly saturates at a level that is almost independent 
of its initial value. 

The V • u terms in the production term for k (see 
Sec. 13. 1|) means that k is also generated behind shocks, 
as can be seen in 4 specific regions in the snapshot at 
t = 0.1 tec'- behind the incident shock sweeping through 
the ambient medium; behind the reflected shock formed as 
the incident shock converges on the axis behind the cloud; 



behind the bow shock formed upstream of the cloud; and 
behind the slow shock driven into the cloud. 

The reflected shock on the axis behind the cloud be- 
comes increasingly oblique as the point of convergence moves 
away from the rear of the cloud, and interacts with the in- 
cident shock to create a do uble Mach reflect ed shock that 
propagates along the axis ([Klein et al.l 1 19941 ) . A powerful 
supersonic vortex ring forms just behind the Mach reflected 
shock, in which a region of high turbulence is generated (see 
the second from left panel of Fig. [H]). While the turbulence 
generated behind shocks decays very rapidly, the turbulence 
associated with the vortex ring is much more persistent, as 
is the turbulence generated at the slip surface around the 
cloud. 

Fig. [8] shows that at later times the turbulence gen- 
erated at the slip surface proceeds to develop into a highly 
turbulent wake with a radius comparable to the initial cloud 
radius. The setup time for the wake is ^ tec, and the sub- 
grid turbulent energy of the cloud material grows and then 
dissipates as the cloud is mixed into its surroundings (see 
Figs. [TTk and [22]) . Note that the core of the cloud has very 
little turbulence associated with it. 

The finite timescale for the development of significant 
turbulence means that simulations with the k-e subgrid tur- 
bulence model with a low initial level of postshock turbu- 
lence produce similar morphologies to those obtained from 
inviscid calculations at early times (t < 0.5 tec)- However, 
the increase in the transport coefficients in regions of high 
turbulence leads to increasing divergence from inviscid cal- 
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culations at later times, and ultimately to a faster destruc- 
tion of the cloud. 

In contrast, a high level of environmental turbulence 
immediately affects the evolution of the cloud, since the 
transport coefficients around the cloud are also high. The 
two right-most panels of Fig. [8] show the highly turbulent 
post-shock flow engulfing the cloud in model c3hil28. Both 
the Umb of the cloud and the bowshock upstream of the 
cloud become broader and less distinct as the high level of 
turbulence leads to strong diffusion across these boundaries 
(see Fig. (6]). Another major difference compared to model 
c31ol28 is that the turbulence downstream of the cloud at 
t = 0.87 tec is roughly as strong as that in the post-shock 
fiow. Hence the turbulent wake which is seen so clearly in 
the centre panel of Fig. [8] is indistinguishable from the sur- 
roundings in the right panel of Fig. [8] Note also that the 
value of k created downstream of the incident shock (see 
the 2nd from right panel of Fig. [8]) is much smaller than 
the initial post-shock value. This refiects the fact that such 
high levels of post-shock turbulence are not naturally gener- 
ated by a shock sweeping through a perfectly homogeneous 
medium. 

The two rightmost panels in Fig. [13] show the turbulent 
energy per unit mass, /c, on a linear scale at t = 1.81 tec 
for low k-e calculations with x = 10^ and different density 
profiles (models c31ol28 and c31osh64). This highlights the 
fact that the strongest turbulence is generated where the 
shear velocity is high, that the central region of the wake 
immediately behind the cloud has a somewhat lower level 
of turbulence, and that further downstream the turbulence 
has penetrated throughout the wake. The opening angle of 
the turbulent layer in model c31ol28 is estimated as ^ 13°, 
which is in ^ood ag reement with experimental results (see 
ICanto fc Ragalll99lL and references therein) , where an open- 
ing angle of ~ 11° is obtained for a Mach 1.3 fiow past a 
stationary medium (the gas behind a Mach 10 shock has a 
Mach number of 1.31). Future work will examine whether 
this agreement with experiment persists as the Mach num- 
ber is varied. 



4.3 Dependence on cloud density contrast 

A range of density contrasts between the cloud and the 
ambient medium is expected. For instance, x '^ 10^ for 
cold atomic clouds embedded in the warm neutral or pho- 
toionzied medium where T '-^ 10^ K, or for warm clouds em- 
bedded in the coronal gas where T ^ 10^ K. For molecular 
clouds embedded in warm gas, x "^ 10^5 while cold atomic 
clouds embedded in coronal gas have x "^ 10^- 

Fig. [9] shows the destruction of a cloud with a density 
contrast x = 10, while Fig. llOlshows the corresponding case 
for a cloud with x — 10^5 both computed with an inviscid 
code. The colour scaling in both figures is identical to that in 
previous figures for easier comparison. The timescale for the 
cloud material to mix into the ambient fiow scales roughly 
with the cloud crushing tiniescale , tec, in agreement wit h 
previous works (JKlein et al.1 Il994l : iNakamura et al.1 [JOQgI V 
However, the normalized growth timescale for RT and KH 
instabilities decreases with increasing x, as is apparent from 
a comparison of Figs. HI [9land[T0l The normalized drag time, 
tdrag/tcc, increases with x, as does the axial stretching of the 
cloud, c/a (see Table[2]). Figs. [T51IT71 also reveal that the ratio 



of the velocity dispersion in the axial to normal directions, 
6vz/6vr, increases with x- These results are all in agreement 
with earlier works. 

The effect of a highly turbulent post-shock fiow is great- 
est at high X (e.g., compare the values of the mean cloud and 
core velocities in Table [2] for the "high k-e" models against 
the "low /c-e" and inviscid models as a function of x)- This 
is because clouds with a high density contrast survive for a 
considerable time after the initial passage of the shock, and 
thus are subject to considerable "buffeting" by the highly 
turbulent postshock environment, whereas at lower values of 
X, the cloud is destroyed relatively quickly after the initial 
passage of the shock. This is also manifest in the increasing 
disparity with x in the evolution of various global quantities 
from the "hi" models on the one hand, and the "no" and 
"lo" models on the other hand, as shown in Figs. 11511171 

4.4 Dependence on cloud profile 

In Fig. [11] snapshots of the density distribution at different 
times are shown for an inviscid calculation of a Mach 10 
shock hitting a cloud with a shallow density gradient (model 
c3nosh64; see Table [1]). While the resolution is such that re 
is equal to the width of 64 cells on the finest grid, the cloud 
in fact extends to r ^ 2 Tc, so the effective resolution is 
similar to the previous models (see Fig. [1]). 

The interacti on of a shock with a smo oth cloud was pre- 
viously studied bv INakamura et al.l (|2006|) for the case where 
X = 10. In this case the cloud offered little impediment to 
the oncoming shock, with the result that the transmitted 
shock and the intercloud shock had similar mean speeds. 
As a result, the intercloud shock did not converge on the 
2;- axis behind the cloud, and the shock compression from 
the downstream side was weak, leading to a slow lateral ex- 
pansion of the cloud. In contrast, the cloud is a much more 
robust obstacle when x = 10^, and we find that it maintains 
many aspects of the evolution seen in sharper-edged clouds 
(c.f. Figs. l4l and [TT]) . Fig. [T2l compares the density structure 
at t = 2.77 tec for inviscid and low k-e calculations (models 
c3nosh64 and c31osh64) . At later times the material stripped 
from the cloud resembles a single tail-like structure, as was 
also the case for a cloud with sharper edges. 

Nevertheless, the shallower density gradient does lead 
to a milder interaction. This is also manifest as a slower 
growth of turbulence around the cloud (c.f. Figs [8] and [13] 
for k-e models with low initial turbulence), with the turbu- 
lent wake not completely forming until t ^ 2 tec- An exact 
comparison of the respective timescales is complicated by 
the fact that the cloud with the shallow density gradient 
is also larger and more massive, and so the true value of 
tec will be different between the models. In fact, we find 
that the forward and rear shocks driven into the cloud con- 
verge at t ^ 0.8 tec in model c3nol28, and at t ^ 1.0 tec 
in model c3nosh64 (we continue to calculate tec using Eq. [2] 
with X — 10^ and re = 1, despite this equation being ap- 
plicable only to clouds with sharp edges). Since these times 
are not too discrepant, it becomes clear that the growth 
of turbulence around the smoother cloud is indeed s lower , 
in agreement with the statement bv INakamura et al.l (|2006r ) 
that it takes more time to form the slip surface around the 
cloud. The maximum turbulent energy per unit mass (i.e. 
kmax) is also higher in model c31ol28. 
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Figure 9. Snapshots of the density distribution from an inviscid calculation of a Mach 10 adiabatic shock hitting a cloud with a density 
contrast of 10 with respect to the ambient medium and with a density profile specified by pi = 10 (model clnol28). The resolution is 
128 cells per cloud radius. The evolution proceeds left to right with t = 0.0, 0.49, 0.97, 1.93, and 3.85 tec- 




u. 



Figure 10. Same as Fig. [9] but for a density contrast of x 
t = 0.0, 0.46, 0.91, 1.82, and 3.65 tec 



10^ (model c2nol28). The evolution proceeds left to right with 
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Figure 11. Snapshots of the density distribution from an inviscid calculation of a Mach 10 adiabatic shock hitting a cloud with a density 
contrast of 10"^ with respect to the ambient medium and with a shallow density gradient, pi = 1 (model c3nosh64). The formal resolution 
is 64 cells per re, but due to the shallow profile the effective resolution is roughly 128 cells per cloud radius. The evolution proceeds left 
to right with t = 0.08, 0.27, 0.46, 0.85, and 1.81 tec- Note the different spatial scale of the plots compared to Fig.|4] 



4.5 A non-uniform post-shock flow 

To obtain further confirmation of the previous results, we 
have simulated the destruction of a cloud hit by a shock 
with non-uniform post-shock flow. To generate the necessary 
random motions of the post-shock flow, perturbations in the 
density and velocity are mapped onto the flow and allowed to 
evolve. These initial pert ubat ions produce normalized stan- 
dard deviations of 0.18, 0.26, and 0.13 in the density, pres- 
sure and velocity of the post-shock flow. The perturbations 
subsequently decay as the shock sweeps up the smooth inter- 
cloud ambient medium and because of the viscosity inherent 
in any numerical code, and the standard deviations of the 
fluctuations noted above decrease. We use a non-AMR setup 
(i.e. a fine, rather than a coarse, G^ grid), since the perturba- 
tions downstream of the shock would cause a large amount 
of grid refinement, thus erasing the benefits that are usually 
possible with AMR. These two issues (decay of the post- 
shock turbulence, and the inability to effectively use AMR) 
further highlight the benefits of using a sub-grid turbulence 
model. 

In Fig.[T3]we show the initial interaction of a shock with 
grid-scale post-shock turbulence and a cloud with a density 
contrast, x = 10- The bowshock around the cloud is much 
less distinct, and the symmetry of the interaction present in 
the previous figures is broken. The shock shows small-scale 
curvature along its surface and is marginally faster than its 
counterpart with smooth post-shock flow (the effective Mach 
number of the shock is about 5 per cent higher, and reflects 
the extra "turbulent" energy mapped into the initial post- 



shock flow). However, the most important difference is that 
the highly turbulent environment destroys the cloud much 
more rapidly than when the post-shock flow is smooth (c.f. 
Figs. [9l and [T4l and also examine Fig. \TE\ , confirming the 
results obtained with the k-e sub-grid turbulence model. 



4.6 Cloud statistics 

Figs. ll5tfT71 show the evolution of various global quantities 
of the cloud as a function of numerical resolution for the 
inviscid and k-e models with x — 10? 10^ ^^^ 10^- ^^ addi- 
tion, the results from a simulation with ^ = 10^ and a shal- 
low density profile are shown in Fig. 1181 Various numerical 
quantities from these simulations are noted in Table (2] As 
already shown, the solutions that the k-e simulations attain 
depend on the initial level of turbulence in the post-shock 
flow which overruns the cloud. For models with low initial 
turbulence, the evolution of the cloud is similar to that at- 
tained from the inviscid code, particularly at lower values 
of X- However, it is clear from the plots in Figs. [TMTHl that 
when the environment surrounding the cloud is turbulent it- 
self, the mixing of cloud material into the surrounding flow 
proceeds at a much faster pace (as shown, for instance, by 
the more rapid growth in the cloud size (see panels a and 
c) and reduction in the core mass (see panel g)) due to the 
enhanced transport and diffusion coefficients. It is also clear 
that the cloud with the shallow density gradient is typically 
less susceptible to high levels of environmental turbulence 
than clouds with sharper boundaries. 
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Figure 12. Comparison between the inviscid (left panel, model 
c3nosh64) and k-e calculation with low initial postshock turbu- 
lence (right panel, model c31osh64) for a Mach 10 adiabatic shock 



hitting a cloud with x 
{pi = 1) at t = 2.77 tec. 



10"^ and a shallow density gradient 



A detailed discussion of the statistics is presented in 
the following subsections. The results obtained from models 
with X — 10^ s^^e discussed initially in each subsection, and 
then a comparison is made to results obtained with lower 
and higher density contrasts and with a shallower density 
profile. Results from the inviscid model with grid-scale post- 
shock turbulence (model clrt32) are also discussed where 
appropriate. 



4.6.1 Cloud shape 

Fig. [T6k-f) shows the time evolution of the rms cloud and 
core radii, a and c, and their ratios, for simulations with 
X — 10^. The transverse dimensions, adoud and acore, de- 
crease during the shock compression stage, and then increase 
during the expansion stage, adoud reaches a maximum of 
^ 4 re (at t ^ 6 tec) in all 3 models, and stays at roughly this 
level until at least t = 12 tee- Although aeore reaches a simi- 
lar maximum, its value subsquently drops precipitously once 
the mass within the core becomes a small fraction (< 15%) 
of its initial value. 

After the initial compression stage, which lasts until 
t ~ tee, the cloud becomes increasingly elongated in the 
direction of the propagation of the shock. The values of Ceioud 



and Ceioud/ttdoud are still growing at t = 12 tee, when the 
simulations were stopped. At this point the cloud material 
is dispersed over a distance of '-^ 15 re in the axial direction 
in all 3 simulations. However, Fig. [10] shows that even at 
earlier times the material stripped from the cloud is highly 
fragmented and does not resemble a single tail-like structure 
(in contrast, the stripped material better resembles a tail 
when X = 10^ - see Fig. [71). 

The axial radius of the core, Ceore, displays the same 
initial behaviour as Ceioud, though its rate of expansion is 
less rapid. Thereafter Ceore behaves in a similar way to aeore : 
it reaches a maximum in all three simulations by t < 8 tee, 
and then declines as material which was formerly within the 
core is mixed into the surrounding flow to the extent that it 
can no longer be identified as "core" material. The values of 
aeore and Ceore cau show abrupt changes as the cloud frag- 
ments, and are generally (though not always) smaller than 
the corresponding values of aeioud and Ceioud- 

The core remains compressed in the axial direction (i.e. 
Ceore/aeore < 1) in simulatious c2nol28 and c21ol28 until 
t ^ 6 tee, after which Ceore /aeore rapidly increases to val- 
ues much larger than unity. In contrast, in model c2hil28, 
Ceore/aeore < 1 at all timcs, and declines to a minimum value 
of 0.32. The core is always less elongated than the shocked 
cloud (i.e. has a smaller value of c/a) in the c2hil28 simu- 
lation, but this is true only for times prior to ^ 6 tee and 
~ 8 tee in the c21ol28 and c2nol28 simulations, respectively. 

Clouds with X = 10 do not expand as much in the lat- 
eral direction, with the values of aeioud and aeore remaining 
below 2 re as the cloud is destroyed (Fig. [T5h and b). In 
contrast, aeioud continues to increase in the latter stages of 
the interaction in model clrt32, due to the strong diffusion 
in this simulation. The enhanced diffusion from the model 
with grid-scale turbulence probably reflects the larger eddy 
sizes in this simulation, and reveals that the process of cloud 
destruction is sensitive to the properties of any turbulence 
in the surrounding medium. 

Ceioud increases to a value of 10 in model clnol28 at the 
time that the simulation is terminated (t = 28 tee), while 
Ceioud /ttcioud asymptotcs at ~ 5.5 at t ~ 25 tee- In model 
cllol28, Ceioud/aeioud follows fairly closely the behaviour of 
model clnol28 until t ~ 17 tee, after which there is a dra- 
matic plunge as mixing of the material reduces hi below the 
threshold for material to be identified as "cloud material". 
In model clhil28, Ceioud /aeioud reaches a much lower maxi- 
mum of 1.5 at t ~ 12 tee- 

Clouds with X = 10^ show different behaviour again 
(Fig.[T7^-e). In models c3nol28 and c31ol28, aeioud and aeore 
are comparable to the initial cloud radius until a consid- 
erable time into the evolution (t < 6 tee)- aeioud eventu- 
ally reaches values of ~ 4 and 6 re in models c3nol28 and 
c31ol28, respectively. In contrast, the core and the more ex- 
tended cloud grow much more rapidly with time in model 
c3hil28, due to the enhanced transport coefficients in this 
case, and the lateral extent of the cloud eventually exceeds 
10 times the initial cloud radius. However, the faster dis- 
persal of the ablated material into the cloud's surroundings 
does not last, and eventually the growth rate of Ceioud and 
Ceore slows below that from models c3nol28 and c3nol28 as 
efficient mixing with the surrounding flow prevents the most 
dispersed material from being identified as material from the 
original cloud. For this reason, the value of Ceioud eventually 
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Figure 13. Snapshots of the turbulent energy per unit mass, k. The three left-most panels show the evolution of k in model c31osh64 
proceeding left to right with t = 0.08, 0.46 and 1.81 tec- In the two rightmost panels a linear scale is used for k at t = 1.81 tec - model 
c31ol28 is shown in the second from right panel, while the rightmost panel shows again model c31osh64. Note also the different spatial 
scale of the plots compared to Fig. [8] 




Figure 14. Snapshots of the density distribution from an inviscid calculation of a Mach 10 adiabatic shock hitting a cloud with a density 
contrast of 10 with respect to the ambient medium and with grid-scale post-shock turbulence (model clrt32) The evolution proceeds left 
to right with t = 0.0, 0.49, 0.97, 1.93, and 3.85 tec- 
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falls below the values obtained in the inviscid and "low /c-e" 
models for each of the cloud density contrasts considered 
(this change occurs at 4.5 < t/tcc ^ 7). This process of mix- 
ing also limits the maximum value of Ccore obtained in the 
simulations. 

Fig. [T8h-e) shows the evolution of the cloud shape for a 
much shallower initial density profile of the cloud but with 
the same x — 10^ ^^ for Fig. [T71 The behaviour of the cloud 
shape is generally fairly similar to the models of the sharper- 
edged cloud (see Fig. [T7)) . 

Figs. [TMTSl and Table [2] show that the values of 
Ccioud/acioud and Ccore/acore generally increase with x and 
decrease with increasing levels of turbulence in the model. 
Since clouds with a high density contrast survive a long time, 
Ccioud/ttcioud can exceed 20 when x — 10^- 

Figs. [T5k), \T6h ) and \T7h) also show the predicted lat- 
eral expans i on of the cloud determined from Eq. 6.1 of 
iKlein et al.l (| 19941 ). This equation is only appropriate for 
t ^tm (i.e. until the lateral expansion of the cloud reaches a 
maximum). A mismatch at early times is also expected, since 
this expression does not properly account for the initial com- 
pression stage. With these caveats in mind, the predicted 
rate of expansion appears to be slightly too fast compared 
to the model results with ;!^ = 10 and 10^, but slower than 
is obtained from the models with x = 10^- The most strik- 
ing observation is that the lateral expansion of the cloud 
in model c3hil28 is very much faster than the theoretical 
prediction over the period 4 < t/tcc ^ 8, which reveals that 
the highly turbulent surroundings diffuse the cloud material 
at a speed faster than the shocked cloud's internal sound 
speed. 



4-6.2 Cloud mass 

Panel g) in Figs. [TSllTSl shows the time evolution of the core 
mass, which reaches half its initial value at t = tmix- The k-e 
models with low initial turbulence show similar evolution to 
the inviscid models. In contrast, the models with high ini- 
tial turbulence show a much faster decline in the core mass. 
We again note that the results from the sub-grid turbulence 
models are much less sensitive to spatial resolution than 
those from inviscid models (Fig. [2}. This is also revealed 
in the mixing time, tmix, which, for instance, is 8.19,8.50 
and 7.93 tec for models c31o32, c31o64, and c31ol28, but 
is 6.81, 11.47 and 8.69 tec for models c3no32, c3no64, and 
c3nol28 (i.e. the latter 3 models show much greater spread 
in tmix). 

Fig. [19] shows how the cloud material mixes with the 
surrounding gas in models c3nol28, c31ol28, and c3hil28. 
The histograms indicate the fraction of cloud mass over a 
range of density bins of width O.lpmax- For pi = 10, about 
80 per cent of the cloud mass has a density ~ Pmax at the 
start of each simulation. At t = 0.38 tec, the shock driven 
into the cloud has increased the density of the gas it has en- 
countered by a factor of ~ 4, while a substantial part of the 
cloud remains unaffected and at its initial density. Eventu- 
ally the transmitted shock sweeps through the entire cloud, 
and densities more than 10 times higher than the initial 
core density occur as the transmitted shock interacts with 
the shock driven into the back of the cloud. At t ~ tec the 
mass distribution function is approximately flat over a wide 
range in density. As the shocked cloud expands the den- 



sity decreases. The turbulent mixing of cloud material into 
the surrounding flow causes the density of the cloud mate- 
rial to eventually approach the shocked intercloud density 
(p ^ 4 X 10~^pmax)- This evolution proceeds slightly faster 
in model c3hil28, as previously noted. 

The double peaked structure of the mass spectrum at 
t = 1.83 tec in mod els c3nol28 and c31ol28 is consistent with 
previous findings (|Nakamura et al.l |2006|) , with the higher 
density peak at p ~ 2.5pmax representing the main core, 
and the lower density peak at p ^ 0.7/>max representing ma- 
terial stripped from it and mixed with surrounding material. 
While there is little difference between the evolution of the 
mass spectrum in models c3nol28 and c31ol28, in contrast 
there is no sign of a double-peaked mass spectrum in model 
c3hil28. This is further evidence that the destruction and 
mixing of interstellar clouds depends on the level of turbu- 
lence present. 

The rate of mass-loss from the core in models c2nol28, 
c21ol28, and c2hil28 can be compared to t he analytical for- 
mula for hydrodynamic ablation given by iHartauist et al.l 
(| 19861 ). To make such a comparison our simulations must 
be appropriately scaled. We therefore assume that the 
cloud is ionized, has a radius re = 2 pc, a core density 
Pc = 4 X 10~^^ gcm~^, and a temperature of 8000 K, 
and is in pressure equilibrium with surroundings of den- 
sity 4 X 10"^'^ g cm~^ and temperature 8 x 10^ K. Both 
the cloud and its surroundings are assumed to have solar 
abundances. The Mach 10 shock then travels at a speed 
of 1360 kms~^ through the ambient medium, heating the 
medium to T = 2.6 x 10^ K, and giving it a velocity and 
Mach number of ~ 1000 kms~^ and 1.3, respectively. The 
rate of mass- loss from a cloud of mass Me through hydro- 
dynamic ablation is Mab = Mc/t ^ /^pit^exp, where t is a 
characteristic destruction/mixing timescale, pi is the char- 
acteristic density of ablated material within a mixing region 
of size l^ around the cloud, and Vexp is the expansion speed 
of material off the surface of the cloud. Momentum conser- 
vation requires that piVe^p ~ pv where p and v are the den- 
sity and velocity of the surrounding flow. If the surrounding 
flow is supersonic, Vexp ~ Cc, since material cannot leave 
the cloud much faster than the sound speed of the cloud, 
Cc. Hence Mab ~ (MeCc)^/^(pv)^/^. With the above param- 
eters, Mab ~ 1-2 X 10~^ M0yr~^ , and the cloud survives 
for approximately 1.6 x 10^ yrs. In comparison the cloud 
crushing timescale, tec = 1.4 x 10^ yrs, so that the cloud 
survives for about 10 cloud crushing timescales before being 
destroyed. 

Fig. I20fe i) compares Mab with the numerically deter- 
mined mass-loss rates, where it is apparent that there is good 
basic agreement between the analytical and numerical mass- 
loss rates. Not surprisingly, the mass-loss rates from the 
numerical models are time-dependent: they decrease as the 
cloud is compressed by the shocks driven into it, and then 
increase significantly as the cloud re-expands. The higher 
mass-loss rates compared to Mab during this latter period 
reflect the increased surface area of the cloud which is not 
taken into account in the analytical theory. The mass-loss 
rate in model c2hil28 is higher than in models c2nol28 and 
c21ol28 for the first t ~ 55, 000 yrs of the interaction - this 
is due to the increased transport coefficients which cause 
the stripping of material at a faster rate. During this time, 
the mass-loss rate from model c2hil28 is also remarkably 
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Figure 15. Time evolution of the cloud shape, mass, volume, mean density, mean axial velocity, velocity dispersion and circulation for 
invsicid and k-e models with a density contrast x = 10, a steep density gradient (pi = 10), and 128 cells per cloud radius. 



constant. The maximum mass- loss rate in all 3 models is 
^ 6x 10~^ Mq yr~^ , or 5 t imes the rate predicted from the 
formula in lHartauist et al.l (| 19861 ). In model c2hil28 this oc- 
curs at the end of the cloud's life. In contrast, the peak rate 
of mass- loss occurs at about 70, 000 yr (t ~ 5 tec), in models 



c2nol28 and c21ol28, which is similar to the time at which 
half of the core mass has been mixed (see Table [2]). In models 
c21ol28 and (especially) c2nol28, the cloud enjoys a more 
sedate ending of its life, with the mass-loss rate declining as 
the core mass decreases. The time at which half of the core 



The turbulent destruction of clouds 21 




t (tec) 



d) ' ' 'J 

c2no128 A 

c2lo128 / \ 




,y .'' '*!.. 

1 1 1 ■ 1 







t (tec) 



c2no128 ■ 
c2lo128 
c2hi128 ■ 





c2no128 ■ 
X c2lo128 ■ 

.'' •-. c2hi128 ■ 



2 4 



t(tcc) 




c2no128 ■ 
c2lo128 ■ 
C2hi128 ■ 





t(tcc) 



Figure 16. As Fig.[T5]but for x = 10^- 



is mixed (tmix) and the destruction time of the core (when 
mcore = 0) are more widely separated in these simulations. 

Fig. 120b) shows an identical analysis for a cloud with 
a smoother density profile (models c2nosh64, c21osh64, and 
c2hish64). The mass of the cloud is 2.9 times greater, and 
the predicted mass-loss rate is twice as high. The cloud is 



therefore expected to survive for 40 per cent longer (i.e. 
230,000 yrs). The initial rate of mass- loss from the sim- 
ulations is approximately 3 times as high as the predicted 
value, which reflects the ease at which the tenuous outer lay- 
ers of the cloud are removed. The variations in the rate of 
mass- loss also appear to be reduced, again consistent with 
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Figure 17. As Fig.[T5]but for x = 10^- 



a milder interaction. Otherwise, similar behaviour to the 
previous case of a sharp-edged cloud is seen, though the 
cloud survives appreciably longer than predicted in models 
c2nosh64 and c21osh64. 

To conclude this section, we note that the cloud in 
model c2hil28 is destroyed in about 65 per cent of the pre- 



dicted time, while in models c2nol28 and c21ol28 the de- 
struction time is in good agreement with the analytical pre- 
diction. Clouds with a smoother density profile survive for 
about 30 per cent longer than the predicted time, though 
if they are overrun by a highly turbulent environment they 
may survive only 60 per cent of the predicted time. A sim- 
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Figure 18. As Fig. [T5]but for x = 10^ cind a shallow density gradient {pi = 1). 



ilar analysis for the simulations with x = 10^ reveals that 
in model c3hil28 the cloud is destroyed more than twice as 
fast as predicted by the analytical theory appropriate for 
a smooth post-shock flow. In model c31ol28 it is destroyed 
in about 75 per cent of the predicted time, while in model 
c3nol28 the destruction time is in good agreement with the 



analytical prediction. The Mach number dependence of the 
mass- loss rate is examined in a subsequent paper (Pittard, 
in preparation). 
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Figure 19. Time evolution of the mass distribution function for models c3nol28 (panels a-e), c31ol28 (panels f-j), and c3hil28 (panels 
k-o). The histograms indicate the fraction of mass contained within each density bin of width O.lpmax, normalized by the total cloud 
mass, rric. From left to right the mass fractions were computed at t = —0.05, 0.38, 0.87, 1.83 and 7.58 tec- 



4-6.3 Cloud volume and mean density 

In each simulation the volume of the cloud core increases af- 
ter the initial compression stage to a maximum value, then 
decreases as the core material is gradually ablated and mixed 
into the surrounding flow. The time of this reversal is typ- 
ically just prior to tmix- The strong turbulence present in 
models c2hil28 and c3hil28 vigorously rips material off the 
surface of the cloud and advects it downstream, so that the 
/^ = 0.5 isosurface is greatly extended. This leads to a rapid 
increase in the volume of the "core" , Vcore , and a commen- 
surate decrease in the core density, pcore (see Figs. [T6h and 
j and I17h and j). In contrast, the destruction of the cloud 
when X = 10 is so rapid that there is not enough time for the 
highly turbulent surroundings to produce much of an effect 
(see Fig. I15h and j), whereas a smooth density profile tends 
to delay the time at which the effects from highly turbulent 
surroundings become obvious (see Fig. [T8h and j). 



4.6.4 Cloud velocity 

The acceleration of the cloud occurs in two stages. The cloud 
is first accelerated to a velocity Vs by the shock driven into it. 
The flow of shocked intercloud gas then further accelerates 
the cloud until they have the same velocity (i.e. 0.725vb 
when M — 10). Since Vs/vh — 1/x^ 5 the second stage 
dominates when x ^ '^•5- If the cross- section of the cloud 
were to remain constant, the characteristic dr ag time for 



sharp-edged clouds would be tdrag,o 



X 



1/2 



(|Klein et al.1 



snarp - 

1 19941 ). However, since the cloud ex pands laterally, th e actual 
drag time is considerably shorter ([Klein et al.lll994l , see also 
Table [2]). 

Fig. 116k ) shows that the acceleration of the cloud is vir- 
tually identical in models c2nol28 and c21ol28 for t ^ 4 tec- 
This is also the case for the core acceleration (Fig. 1161). In 
contrast, the stronger turbulent mixing in model c2hil28 in- 
creases the momentum coupling of the surface of the cloud 
with the ambient flow so that the cloud accelerates more 



quickly to the post shock intercloud speed. This effect is even 
more pronounced in model c3hil28 (see Fig.fTTk). The accel- 
eration of clouds with smooth density profiles is more gentle 
(cf. Figs. [TTk and [T8k) . as expected. In all cases the cloud 
accelerates faster than the core (i.e. (vz, cloud) > ('^z,core)), as 
expected. 

iKlein et al.1 (|l994l ) and iNakamura et al.1 (I2OO6I ) pre- 
sented a semi-analytic theory for the acceleration of sharp- 
edged an d smooth- edged clo uds, respectively. The formula- 
tion from iKlein et al.l (| 19941 ) converted to the rest frame of 
the initial cloud is shown in Figs. [T6k ) and [T7k) . where the 
value of tm obtained from models c2nol28 and c3nol28 is 
used. In each case, the results obtained with three different 
valu es for the dr a g coe fficien t are shown: Cp = 0.2, .5, and 
1.0. iKlein et al.1 (|l994l ) and INakamura et al.1 (I2OO6I ) adopt 
Cd = 1.0, but it is not clear that this is the best choice. For 
solid spheres, experime nts have determined tha t Cd = 1.0 
occurs when Re ~ 100 ([Landau fc Lifshitj|l959l ), which is a 
far lower Reynolds number than occurs in the astrophysical 
settings that we are interested in (see Section [2?2)) . However, 
interstellar clouds are compressible and are ablated in shock- 
cloud interactions, so one would not necessarily expect the 
drag coefficient to be similar to the solid sphere case. Never- 
theless, clouds with a high density contrast should, at least 
initially, behave somewhat like solid spheres. If this is the 
case, more appropriate values for Cd would be ^ 0.5 when 
2 X 10^ < Re < 2 X 10^ and ^ 0.2 for Re > 2 x 10^ 

In fact, we find that there is poor agreement between the 
theoretical acceleration and the results from model c2nol28 
for values of Cd ^ 1 (Fig. flBk). The closest match is ob- 
tained with Cd ~ 1.0, though the theoretically estimated ac- 
celeration is too high for t < 4 tec, and too low at later times. 
This discrepancy is the opposite of that noted in lKlein et al.l 
(1994). The poor match is at least in part due to the poor 
agreement between the predicted and model evolution of 
ttcioud (Fig. flBh). Since Cd ~ 0.2 — 0.5 for high Reynolds 
number flow over a solid sphere, we conclude that clouds 
with X = 10^ are a poor approximation to a hard sphere. 
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Figure 20. a) Mass loss rate for models c2nol28, c21ol28, and c2hil28 compared to the mass-los s formula of iHartauist et al.l (|l986h . b) 
Mass loss rate for models c2nosh64, c21osh64, and c2hish64 compared to the mass-loss formula of lHartquist et ahl (|l986l ). In both panels 
t = corresponds to the time at which the intercloud shock is level with the centre of the cloud. 



In contrast, Fig. [T7k ) shows that the theoretical accel- 
eration of the cloud agrees much better with the numerical 
results from models c3nol28 and c31ol28 where x — 10"^ • 
The best agreement is obtained with Cd ~ 0.5, which is con- 
sistent with the assumption of fully developed turbulence in 
the k-e turbulence model. We therefore conclude that clouds 
with X = 10^ are a better approximation to a hard sphere. 



4-6.5 Cloud velocity dispersion 

The interaction of shocks with clouds produces substan- 
tial vorticity and velocity dispersion, which may be a key 
mech anism for g enerating turbulent moti ons in the ISM 
(e.g., iKornreichfc Scalo 2000; Mac Low fc KlessenI IJOOl ) . 
Although the sub-grid turbulence model deals with turbu- 
lent motions and mixing on scales smaller than the cell size 
of the numerical grid, larger scale turbulent motions can be 
directly measured from the velocity dispersions in the ax- 
ial and radial directions, 6vz and 6vr respectively. The time 
evolution of these quantities is shown in panels m) and n) 
of Figs. [TMTSl 

Our results for model c2nol 28 can be compared to 
Figs. 9 and 10 in lKlein et all ([1994'). The level of qualitative 
and quantitative agreement is basically good. In particular, 
we confirm that the radial velocity dispersion is generally 
somewhat less than the axial velocity dispersion. However, 
we find that the maximum values of Svz and Svr occur some- 
what later in our simulations. Comparing models c2nol28, 
c21ol28, and c2hil28, we see that 6vz is more sensitive to a 
higher initial level of turbulence in the environment around 
the cloud than 6vr. 

A comparison with the x = 10 and 10"^ simulations re- 
veals that Svz/Svr increases w ith x, as is also apparent in 
the work of iKlein et al.l (|l994r ). 6vz peaks at about the same 
value in c2hil28 and c3hil28, indicating that the random 
instability-induced grid-scale motions which develop are lim- 
ited by the high viscosity introduced by the subgrid turbu- 
lence model in these simulations. 6vz peaks at higher values 
in models c3nol28 and c31ol28, compared to models c2nol28 
and c21ol28, due to greater growth of RT and KH instabil- 
ities resulting from the longer drag and mixing time in the 



former models. The velocity dispersion appears to be slightly 
reduced for clouds with smooth boundaries, with some sign 
that the time evolution of Svz attains a broader maximum. 
In all simulations, dvr peaks slightly prior to the time of the 
maximum radial extent of the cloud, tm- 

Fig. [21] shows the mass distribution function in veloc- 
ity space within the entire cloud integrated along the z- 
axis for models with x = 10"^. The histograms indicate the 
mass contained within a corresponding velocity bin with a 
width of 0.01 i^ics, where Uics is the postshock ambient ve- 
locity. The shock initially driven into the cloud has a speed 
^h/\/x = 0.042 Uics- By t = 3.74 tec the majority of the 
mass within the cloud in models c3nol28 and c31ol28 has 
been further accelerated by additional shocks and momen- 
tum transfer from the surrounding flow to roughly twice this 
speed. In contrast, the material stripped off the surface of 
the cloud moves at speeds up to Uics, with a small frac- 
tion exceeding Uics as a result of turbulent motions. Further 
stripping and momentum transfer results in a gradual accel- 
eration of the majority of the cloud material to speeds near 
i^ics as time progresses. Not unexpectedly, the acceleration 
of cloud material in model c31ol28 is more rapid than in 
model c3nol28, and even more so in model c3hil28, where 
even the slowest moving cloud material has a speed in excess 
of 0.6 i^ics at t ^ 12 tec. 



4-6-6 Cloud vorticity 

A key feature of the interaction of a shock with a cloud is the 
development of powerful vortex rings. Vorticity (cj = V x v) 
can be produced (or destroyed) when pressure and den- 
sity gradients are not aligned (i.e. by a curl in the accel- 
eration), and by viscosity. The c orresponding circulation, 
r = Jcj • dA. iKlein et aJ] (| 19941 ) showed that the vortic- 
ity production can be classified into 4 components. Vortic- 
ity at the interface between the cloud and the surrounding 
flow is produced by the initial passage of the shock (Pshock) 
and the subsequent postshock flow (Ppost). The third com- 
ponent is connected with the triple points associated with 
the Mach-reflected shocks behind the cloud (Fring), while 
the fourth component is the vorticity produced in the cloud 
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Figure 21. Mass distribution as a function of Vz for models c3nol28 (panels a-c), c31ol28 (panels d-f), and c31iil28 (panels g-i) at 
t = 3.74, 7.58 and 11.9 tec- The histograms denote the mass contained within a velocity bin of width 0.01 Uics for a line of sight parallel 
to the z-axis. The integrated mass is the mass of the cloud, rric. 



(rcioud) a nd is smaller than the ot her components by a factor 
^ X~^^^ • iiNakamura et alJ (|2006h showed that the vorticity 
produced by a shock overunning a cloud with smooth bound- 
aries is qualitatively si milar to that prod uced from clouds 
with sharp boundaries. iKlein et al.l (| 19941 ) constructed ana- 
lytical expressions for 3 of the above components which we 
reproduce below: 



(26) 

I'brc, (27) 

(28) 



Tpost is the dominant c ontribution to the circulation when 
X is large. iKlein et aU (|1994) did not attempt to model 
Tcioud, since this component is small. For clouds with smooth 
boundaries, Fpost must be further multiplied by (1 — x~^)^ 
(JNakamura et al.ll2006l l. 

The time evolution of the circulation from our numer- 
ical models is shown in panel o) of Figs. 11511181 Table [3] 
lists values estimated for Fpost and the total circulation, 
Ftot = Fshock + Fpost + Fring • lu cach case we use the values 
for tdrag in Table [2] Obviously, the circulation computed us- 
ing the cloud and core drag times is similar when tdrag for 
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the cloud and the core are also similar (see, e.g., models of 
sharp-edged clouds with ^ = 10 and 100). However, when 
tdrag for the cloud and the core differ (as occurs for mod- 
els with a shallow cloud density profile), we find that the 
numerically determined circulation is in much better agree- 
ment with theoretical estimates if the drag time for the core 
is used. Figs. fTMTSl indeed show that the mean velocity of 
the cloud and core as a function of time become increasingly 
disparate with increasing x- Clearly it is the motion of the 
core with respect to the postshock ambient medium which 
dominates the generation of circulation - after all, this is 
where the highest velocity shear occurs. 

Table [3] reveals that the model results and the theoret- 
ical estimates are generally in good agreement. The contri- 
bution of the 3 main components is best seen in Fig. ll5bV 
The initial rise to maximum is caused by the initial passage 
of the shock (Fshock), while the subsequent drop is caused 
by the formation of the supersonic vortex ring behind the 
cloud, Fring- The increase in circulation after this minimum 
(until t ^ 2 tec) reveals the vorticity production due to the 
post-shock flow over the cloud. 

The total circulation in the models is often in excellent 
agreement with the predictions (for instance, the peak cir- 
culation in models c31ol28, c3hil28, c31osh64, and c3hish64 
are 33.8, 23.3, 65.2, and 40.0, while the theoretical predic- 
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Table 3. Theoretical estimates for the total circulation and the 
component produced by the postshock flow. In each case the cal- 
culated value uses the drag-time for the cloud (core), as given in 
Table E] -Fshock = 1.54,2.03, and 2.18 when x = 10, 10^ and 
10^, respectively, while the value of Fring is always 0.75. 



Model 


1 post 


— Ttot 


clnol28 


0.45 


(0.49) 


1.24 


(1.28) 


cllol28 


0.47 


(0.49) 


1.26 


(1.28) 


clhil28 


0.32 


(0.38) 


1.11 


(1.17) 


c2nol28 


4.30 


(4.35) 


5.58 


(5.63) 


c21ol28 


4.33 


(4.35) 


5.61 


(5.63) 


c2hil28 


3.47 


(3.73) 


4.75 


(5.01) 


c3nol28 


29.2 


(42.1) 


30.6 


(43.5) 


c31ol28 


30.4 


(31.8) 


31.8 


(33.2) 


c3hil28 


20.4 


(21.8) 


21.8 


(23.2) 


c2nosh64 


4.61 


(6.28) 


5.89 


(7.56) 


c21osh64 


4.56 


(4.93) 


5.84 


(6.21) 


c2hish64 


3.51 


(4.39) 


3.08 


(3.85) 


c3nosh64 


28.0 


(41.0) 


29.4 


(42.4) 


c31osh64 


31.5 


(60.2) 


32.9 


(61.6) 


c3hish64 


27.3 


(34.9) 


28.7 


(36.3) 



in time, while in model c3hil28 the kinetic energy exceeds 
the thermal energy during the period 4 < t/tcc ^ H- In 
model c3hish64, the kinetic energy always exceeds the ther- 
mal energy. It is clear that the cloud material has still not 
fully mixed into the surrounding flow at the end of most of 
the simulations, since the thermal energy often significantly 
exceeds the kinetic energy at this point (obvious exceptions 
are models c2hil28, c3hil28, and c3hish64). 

The fraction of the cloud energy in small-scale subgrid 
turbulent motions increases with x for the "low /c-e" models, 
with maxima of 0.0062, 0.017, and 0.029 for models cllol28, 
c21ol28, and c31ol28, respectively. When the environment 
is highly turbulent, the turbulent energy fraction is 0.079, 
0.107 and 0.157 for models clhil28, c2hil28, and c3hil28, 
respectively. The latter value is similar to the fractional en- 
ergy in turbulence in the post-shock flow prior to its im- 
pact with the cloud, and demonstrates how high levels of 
upstream turbulence may be maintained through sequen- 
tial shock-cloud interactions. Clouds with smoother density 
profiles produce a slightly smaller peak fractional energy in 
turbulent motions (0.020 and 0.140 for models c31osh64 and 
c3hish64, respectively). 



tions are 33.2,23.2,61.6, and 36.3, respectively). In some 
models (e.g., c2nol28, c21ol28, c3nol28 and c3nosh64) the 
circulation is somewhat higher than predicted. Spurious vor- 
ticity can be g:ene rated at the boundary of grids of different 
refinement levels (|Plewa fc MiillerlboOlh , but this does not 
seem to be the case here (since, for instance, one would ex- 
pect a larger discrepancy for model c31ol28 than for model 
c21ol28). 

The total circulation in models with x — 10^ tends 
to reach a lower maximum when there is high postshock 
turbulence, due to the more rapid destruction of the cloud 
in such circumstances. Models where the cloud has a shallow 
density profile generate more circulation, principally because 
of the larger cloud mass and the longer drag times in these 
models. The decay in Ftot at later times in some models may 
reflect the inherent viscosity diffusing vorticity as turbulence 
is dissipated into heat. 



4-6.7 Energy fractions 

As the cloud material accelerates it gains kinetic energy, 
while its thermal energy increases due to shock heating, 
adiabatic compression, and heat transport from turbulent 
mixing with the hotter surrounding flow. There is also likely 
to be some numerical transport of heat. The turbulent en- 
ergy increases as shear motions generate a turbulent cascade. 
Eventually, the cloud material should acquire the same ki- 
netic and thermal energy density as the ambient medium, 
and the turbulent energy should dissipate as heat. For a 
strong adiabatic shock, the ratio of kinetic to thermal en- 
ergy in the postshock flow, E]^/Eth = 7Mpg/3, where Mps is 
the postshock Mach number. In our simulations Mps = 1.32, 
so at late times we expect that Ei^/Eth — 1. 

Fig. [22] shows the time evolution of these energies and 
the fraction of the total energy in small-scale (subgrid) tur- 
bulent motions for simulations with x — 10? 10^ and 10^. In 
models c3nol28 and c31ol28, the thermal energy of the cloud 
generally exceeds its kinetic energy at any particular instant 



5 DISCUSSION 

The interaction of shocks, winds and jets with clouds, and 
the collapse and/or destruction of the cloud, is of fundamen- 
tal importance to studies on star formation, the ISM, feed- 
back and galaxy evolution, and the evolution of diffuse as- 
trophysical sources, such as planetary nebulae, wind-blown- 
bubbles, supernova remnants, Hii regions, galactic winds, 
and the intracluster medium. 

This work has examined the interaction of a shock with 
a cloud. However, clouds with a high density contrast easily 
survive the initial passage of the shock, and then find them- 
selves immersed in a subsonic or mildly supersonic post- 
shock flow which to all intents and purposes resembles a 
wind of the same Mach number. Hence the simulations pre- 
sented in this work are also relevant to scenarios where 
clouds are ablated by a wind. In a future paper we will 
examine in detail the differences between shock-cloud and 
wind-cloud interactions. 

In the following subsections we discuss shock-cloud in- 
teractions in supernova remnants (SNRs), and the tails 
formed behind clouds in various types of wind-cloud inter- 
actions. Jet-cloud interactions, though not discussed in this 
work, can be seen in Herbig-Haro objects and AGN. In some 
sources, clouds are more commonly referred to as bullets, 
clumps, globules, or knots. 



5.1 SNR-cloud interactions 

The study of shock cloud interactions has typically focussed 
on SNRs. In most cases the shocks driven into the clouds are 
radiative, so such interactions are not directly comparable 
to the adiabatic simulations presented in this paper. Nev- 
ertheless, in recent years there has been much progress in 
unravelling the nature of these interactions and it is worth 
discussing some of their interesting features. 

An interesting observation is that some clouds show 
clear signs of a bowshock (for instance, the Cynus Loop's 
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Figure 22. Time evolution of the cloud kinetic and thermal energies for invsicid and k-e models. The cloud turbulent energy at subgrid 
scales and its fraction of the total cloud energy is also shown for the k-e models. 



southeast clou d (IGraham et al.lll995|l, and also cloud s on the 
western limb ([Levenson, Graham h WaltersI I2OO2I I), while 
others do not (e.^., the XA clou d in the Cygnus Loop 
(jDanforth, Blair fc Raymondll2001), the sou thwestern cloud 
of the Cygnus Loo p ([Patnaude et al.ll2002i ), and Vela FilD 
(JMiceli et al.ll2006r )V Possible explanations for t he lack of a 
bows hock are the earliness of the interaction (^Klein et alJ 
199J), a smooth cloud density profile (Nakamura et al^, 
20061 ). or a high ellipticity of the cloud iMiceh et al.ll2006l ). 



Therm al conduction can a lso reduc e the visibility of bow 
shocks (jOrlando et al.ll2005l ) (though [Marcolini et all (l2005l ) 
show that in the case of a wind-cloud interaction it can en- 
hance its visibility), while co smic rays can smooth shocks 
out (e.g., IWagner et al.ll2006l ). In addition to these mecha- 
nisms, we note that a highly turbulent post shock flow may 
also hinder the formation of a clearly defined bowshock (see, 
e.g.. Fig. I14p . Variations in pre-shock density due to inter- 
stellar turbulence have been estimated to be ^ 20 per cent 
([Raymond et al.l boOTJ ). For SNRs, post-shock turbulence 
could be generated by pre-shock interstellar turbulence, or, 
at least in young SNRs like SN 1006, by clumps of ejecta 
and/or long RT fingers affecting the contact discontinuity 
and blast shock. Parti cle acceleration also p roduces turbu- 
lent motions (see, e.g.. |jones h Ellisonlll99lL and references 
therein) . Finally, we note that although thermal conduction 
has been invoked as necessary to explain the hot corona 



surounding the Vela FilD cloud (JMiceh et al.ll2006l ). it is 
possible that the turbulent transport of heat could instead 
account for this. 

It is also interesting to note that all of the numerical 
simulations of shock-cloud interactions have assumed that 
the interaction takes place in the so-called "small-cloud" 
limit, where the size of the cloud is sufficiently small that 
the properties of the post-shock gas do not change signifi- 
cantly in the time that it takes for the cloud to be crushed or 
destroyed. This is unlikely to be the case for denser clouds. 
In addition, the interaction of a superbubble with a cloud is 
better described as a shell-cloud interaction, and will be the 
subject of a forthcoming work. 



5.2 Wind-cloud interactions: tails 

Material stripped off clouds by the passage of a shock 
or fast wind is frequently manifest as a long tail behind 
the cloud (see, e.g.. Fig. [7]). Although this work does not 
specifically focus on tails, the tails themselves are inter- 
esting as their properties may allow a diagnosis of the 
flow past the clo ud (e.g., Dvs on. Hartquist fc Bird Il993l : 
IPvson et al.ll2006h . Examples of resolved tails include those 
in NGC 7293 (the Hel i x nebula; lo 'Dell Hennev h Ferlandl 
I2OO5I : iHora et al.ll2006l : iMatsuureTet al.ll2007l ). the complex 
sub-structure of knots and wakes seen in the Orion Molecu- 
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lar C l oud OMCl (e.g., [Allen fc Burtonlll993l: ISchultz et alj 



19991: iTedds Brand fc Burtonl Il999l : iLee fc BurtonI l2QQ0l : 



Kaifu et al 



2OOOI ) 



and the radial protrusions from 
the bipolar nebula surrounding the super massive star 
n Carinae (IWeis. Duschl fc Chul il999l : ICurrie et all I2OOOI : 
iRedman. Meaburn Sz Hollowavl I2OO2I ) . Comet-like tails are 
also seen extending from the Galactic Centre source 
IRS 7, a red supergiant (lYusef-Zadeh fc MorrisI Il99ll : 



4 



ISerabvn. Lacv Sz Achtermannlll99l | ). and from behin d Mira, 
an asymptotic giant branch star ([Martin et al.ll2007l ). 

Tail- like structures are also see n in galactic winds 
([Cecil et al.ll200ll : lOhvama et alJ 12002). and resemble simi- 
lar features seen in simu lations (e.g., IStrickland Sz Stevens! 
l200d : ICooper et all l2008h . Observations and simulations of 
galactic winds suggest that the clumps at the heads of the 
filaments are material which is ripped out of the galactic disk 
as the galactic wind develops. Beautiful filamentary struc- 
tures are also seen in some galaxy clusters, o f which the most 
prom inent example is in the Perseus cluster ([Conselice et al.l 
|200l|). Numerical simulations indicate that a high veloc- 
ity wind is necessary for these to form, otherwise mate- 
rial stripped from the clou d sinks under grav ity almost as 
quickly as the cloud itself ([Pope et al.l 120081 ). Tails in the 
intracluster medium may help to regulate the heating of 
the central regions of galaxy clusters by dissipating some 
of the energy injected by the central active galactic nucleus 
(AGN), and therefore better couple the AGN to the intra- 
cluster medium. 

The Reynolds number is likely to be high enough in all 
of the above cases that the tail is fully turbulent. Indeed, 
the emission in Mira's tail can be explained if molecular hy- 
drogen is excited by the t urbulent mixing of c ool molecular 
gas and shock- heated gas ([Martin et al.ll2007l V 



6 SUMMARY AND CONCLUSIONS 

This is the first in a series of papers investigating the tur- 
bulent destruction of clouds. Here we have investigated the 
destruction of a cloud by an adiabatic shock using a hy- 
drodynamical code which incorporates a subgrid k-e turbu- 
lence model, in which an attempt is made to calculate the 
properties of the turbulence and the resulting increase in 
the transport coefficients. The results are compared against 
those from a hydrodynamical code which solves only the in- 
viscid Euler equations of fluid motion. The motivation for 
this study is that fully developed turbulence is prevented 
by the artiflcial viscosity in all numerical codes, though it is 
expected in astrophysical environments where the Reynolds 
number of the shock-cloud interaction is ^ 10^. The effect 
of a highly turbulent environment sweeping over the cloud is 
also investigated - all other works studying the interaction of 
a shock with a single cloud have assumed that the post-shock 
environment is completely smooth, though this is clearly not 
the case when there are upstream inhomogeneties. Our main 
results are summarized thus: 

(i) The evolution and destruction of clouds with inviscid 
and k-e models occurs at roughly the same speed when the 
post-shock flow is smooth. This is because the turbulence 
takes some time to be generated and for its effects to be 
subsequently felt. However, it is clear that there are increas- 
ing differences between such models as the density contrast. 



X, between the cloud and inter-cloud medium increases. We 
also show that the k-e model results are far less dependent 
on the spatial resolution than calculations with an inviscid 
code. This behaviour is attractive given the recent interest 
in simulations with multiple clouds. 

(ii) Turbulence is mainly generated at the slip surface 
around the cloud, though a small amount is also generated 
behind shocks. The set-up time to form a turbulent wake 
behind the cloud is ^ tec- There is slower growth of tur- 
bulence around clouds with a smooth density proflle. The 
opening angle of the turbulent mixing layer is consistent 
with experimental results. The fraction of energy in small- 
scale sub-grid turbulent motions initially increases with time 
as turbulence continues to be generated by the velocity shear 
of flow past the cloud, but then decays as turbulent energy 
subsequently dissipates as heat. The peak value of the tur- 
bulent energy fraction increases with x due to the longer 
drag-time of denser clouds and the higher velocity shear at 
the cloud surface, and decreases as the boundary of the cloud 
becomes smoother. The interaction of shocks with smooth 
clouds is generally milder, as previously reported. 

(iii) Clouds which are subject to a highly turbulent post- 
shock environment are destroyed significantly quicker than 
those within a smooth flow, due to the enhanced transport 
and diffusion coefficients. This effect increases with x, since 
high density clouds are subject to a longer period of "buf- 
feting". Clouds with small density contrasts (e.g., x — 10) 
are destroyed so quickly after the initial shock passage that 
they experience very little "buffeting" . Strong environmen- 
tal turbulence increases the momentum coupling between 
the cloud and its surroundings, which increases the acceler- 
ation of the cloud. The degree by which the destruction of 
the cloud speeds up should depend on the strength of the 
turbulence imposed on the post-shock flow, but here our in- 
tention is simply to demonstrate that it does, in fact, occur 
more rapidly. 

(iv) The rate of mass-loss from a cloud overrun by a 
Mach 10 shock due to hydrodynamic ablation is found to be 
broadly consistent with theoretical expectations, but shows 
large variations over the destruction period: the peak mass- 
loss rate is about 5 times higher than the time-independent 
value from a simplifled theory. Clouds in a highly turbulent 
environment can be destroyed approximately twice as fast. 
However, signiflcant differences between theoretical and nu- 
merically determined mass-loss rates exist at higher Mach 
numbers (Pittard, in preparation). 

(v) Material stripped off the cloud is fragmented and ir- 
regular in structure when x ^ 10^- In contrast, the mass- loss 
from models with x — 10^ better resembles a single tail-like 
feature. The length-to- width ratio of the tail increases with 
X (denser clouds live longer, allowing material to be more 
dispersed in the axial direction) , but decreases with the level 
of environmental turbulence. 

(vi) The vorticity generated in the interaction is broadly 
similar in models with differing levels of post-shock turbu- 
lence when X ^ 10^- However, models with high environ- 
mental turbulence generate substantially less total circula- 
tion when X = 10^ compared to models with a smooth post- 
shock environment, due to the more rapid destruction of the 
cloud. This may self- limit the total circulation that can be 
generated in such interactions. 

(vii) Conflrmation of the general speeding up of cloud 
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destruction in a highly turbulent post-shock environment is 
attained by additional calculations where strong, grid-scale 
motions and inhomogeneities are imposed in the post-shock 
flow of an inviscid calculation. Differences in the resulting 
evolution compared to those from the subgrid turbulence 
models with high post-shock turbulence reveal that the in- 
teraction is sensitive to the details of the turbulence, such as 
the maximum eddy size. The degree by which the destruc- 
tion speeds up should again depend on the strength and 
properties of the turbulence imposed on the post-shock flow. 
Future simulations with a non-smooth inter-cloud medium 
would also be of interest. 



7 FUTURE WORK 

Irrespective of the shortcomings of the k-e model (see, e.g., 
lDavidsonll2004r ) , it is clear that in shock-cloud interactions, 
turbulence (both pre-existing and newly generated) plays an 
important role, and adds a new dimension to the parameter 
space that has hitherto been studied. An obvious extension 
of this work is to three dimensions: previous work has shown 
that many features seen in two-di mensional simulation s are 
unstable in three dimensions (St one Sz Normanlll992l V Fu- 
ture work will also examine the Mach number dependence 
of shock-cloud interactions and the interaction of winds and 
dense shells with clouds. 

Synthetic signatures of the interaction should also be 
compared against observations. A comparison of the syn- 
thetic emission from models against X-ray o bservations has, 
for in stance , already been carr ied out by iMarcolini et al.l 
(|2QQ5l ) and iMiceh et all (|20Q6h . Other observational sig- 
natures also need to be examined in greater detail: 
IWestmoauette et al.l (|2007al lb[) have recently concluded that 
the broad emission wings to Ka line profiles which are ob- 
served throughout the central regions of starburst galax- 
ies arise in a turbulent boundary layer at the interface be- 
tween hot ga s flowing past cold gas stripped from clouds 
([Melnick et al. 1999 ) , but the range in velocities may instead 
reflect the acceleration of material along the tail, rather than 
the turbulent motions within the mixing layer. 

In many environments flows interact not with a sin- 
gle cloud, but with many clouds. Hydrodynamic simulations 
of the complex interactions between multiple clouds and a 
tenuo us flow have been s t udied bv Ijun, Jones Sz NormanI 
([1996'), 'Poludnenko et al.' ('2002'), 'Steffen k Lopez' ('200^ 
[PiTt ard et al. (2005), Melioh, de Gouveia Dal Pino & Ra^ 
(EOOSI) . iTenorio-Tade et al.1 (l2006l). [Sutherland fc Bicknell 
J2OO3), lYirak et all (I2OO8I ) and ICooper et al.1 (I2OO8I ). Most 
of these works focussed on the changes to the global flow, 
but the interaction between two or more long-lived c louds in 
close p roximity has been examined in detail bv iPittard et al.l 
(|2005l ). In all situations, mass injection into a flow due to 
the destruction of clouds enhances the thermal pressure of 
the flow at the expense of the flow's ram pressure. If the 
flow is slowed and pressurized enough, it may induce the 
gravitational collapse of clouds and trigger new star forma- 
tion. This process could be a central mechanism for feedback 
in the interstellar medium (e.g., in starburst regions), and 
a self- consistent hydrodynamical model of it is a long-term 
goal. 
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